Data type of ITensor

Question

Hello,

I considering to use ITensor as a part of generic finite element lib (mofem.eng.gla.ac.uk). I like to use ITensor with automatic differentiation provided by lib adol-c.

It is easy way to use data type diffrent than real (double) or complex and use type defended buy adolc, f.e. instead double -> adouble.

Can you point me to revenant documentation.

Kind regards,
Lukasz

Answer 1

Hi Lukasz,
I'm not sure of the best answer to your question. First, you'd have to tell me more about what this other data type is like. I'm not sure what "adolc, f.e" means, for example.

But I can partially answer your question at a couple different levels:

1. The two types for numbers used in the library are Real, which is a typedef (alias) for double and Cplx which is a typedef for std::complex. You could always try changing the definition of Real from double to some other type (in the file itensor/types.h). This would probably leave a large fraction of the code working, but would certainly prevent other parts of the code from compiling. The reason is that there are places where we use library functions (such as std::fabs) which work for types like double, float, etc. but which may not be defined for the type you want to support. But it might be easy to fix these parts of the code to work for your type.

2. Another higher-level approach you could take would be to define a new ITensor storage class based around your data type. The usual storage type of an ITensor is Dense where T=Real or T=Cplx. You could create a new ITensor storage type "ElementStorage" or something and register it with the storage system. This type could be modeled on the Dense type but customized for your needs. I haven't written much documentation on how to do that because it's so new and I have other things to document first. But it's pretty easy to do once you see an example of how.

Probably #2 is your best bet.

Let me know which approach sounds closer to what you're trying to do and we could discuss more about how to proceed.

Comments

Thanks for swift answer.

ADOL-C is a library https://projects.coin-or.org/ADOL-C which is used to do (automatic) differentiation.

In my case I have tensorial field and I would like to calculate some derivatives over that field. For example I have a stress tensor which is function of gradient of defamation (tensor as well) and like to calculate derivative stress over gradient of deformation.

In principle I can calculate derivatives  by hand, but some equations are complex and resulting equations are long. Using ADOL-C (or other type of software), I can calculate derivatives exactly and avoid long derivation and then implementation. It is very convenient for many users.

It looks that second option is the best. The trick how ADOL-C (and other of this type libs) works is that inside of double other type is used f.e. adouble, or activedouble ....

I will have a go with 2nd option, and see if it works.

---

Yes I think the 2nd option could be what you need.

In lieu of a proper documentation of the storage system, here is a brief tour of how it works (using examples from the ITensor code):

* Make your new storage type (it doesn't have to do anything yet) and register it with the storage system: edit the file itensor/itdata/storage_types.h and follow the instructions in that file, then recompile the whole library.

* To make a custom ITensor with your new storage type do:

    auto is = IndexSet(i1,i2,i3); //these are the indices it will have
    auto dat = MyDataType(...); //make the storage type, it can be whatever type
    auto T = ITensor(std::move(is),std::move(dat)); //create the ITensor with custom storage

* To see how the "doTask" system for calling methods on various storage types works, let's use the example of computing the norm and look at the Dense<Real> storage type.

  1. First the norm(T) function is called, code at itensor_interface.ih line 806.

  2. This function calls `doTask(NormNoScale{},T.store());` where T.store() is an "opaque" storage pointer i.e. we don't yet know what storage type the ITensor actually has.

  3. The doTask function "unwraps" the storage pointer and discovers the storage type; this is handled automatically for you.

  4. Assuming the actual storage type is Dense<Real>, the doTask(NormNoScale{},Dense<Real>) function is called, defined on line 125 of itensor/itdata/dense.cc.

  5. This function recieves a Dense<Real> object as a const reference and can then compute the norm of this object as appropriate. doTask functions can return any type you want, as long as all overloads of doTask which take the same "task object" (such as NormNoScale, Contract, PlusEQ, etc.) all return the same type.

For another example of the doTask system at work look at this short article on the ITensor website:
http://itensor.org/docs.cgi?page=formulas/extractdense

Here is an article talking about the concept of the doTask system:
http://itensor.org/docs.cgi?page=articles/storage

One thing that is helpful about the ITensor storage system is that you don't have to define all possible operations on a new storage type for the code to compile and for you to begin using it. You can start by defining a few operations, such as doTask(NormNoScale, MyType), doTask(Contract,MyType) etc. and only use those to test. Then when things are working you can continue to define more operations like doTask(PlusEQ,MyType) etc.