"""
Fortran module Tensors: June 2017 (dj)
Containing basic operations for tensors as
allocation, reading writing.
**Authors**
* D. Jaschke
* M. L. Wall
**Details**
The error status can be decoded the following way
* 1: allocation / deallocation problem
* 2: rank mismatch
* 3: dimension mismatch
* 4: unknown argument
* 5: argument not valid
The following subroutines / functions are defined for the
applicable data type
+----------------------+-------------+---------+
| procedure | include.f90 | mpi.f90 |
+======================+=============+=========+
| addproject | X | |
+----------------------+-------------+---------+
| add_dummylink | X | |
+----------------------+-------------+---------+
| conj | X | |
+----------------------+-------------+---------+
| copy | X | |
+----------------------+-------------+---------+
| create | X | |
+----------------------+-------------+---------+
| create_id | X | |
+----------------------+-------------+---------+
| dagger | X | |
+----------------------+-------------+---------+
| destroy | X | |
+----------------------+-------------+---------+
| dot | X | |
+----------------------+-------------+---------+
| equals | X | |
+----------------------+-------------+---------+
| expm | X | |
+----------------------+-------------+---------+
| expmh | X | |
+----------------------+-------------+---------+
| finalize (splitlink) | X | |
+----------------------+-------------+---------+
| fuse | X | |
+----------------------+-------------+---------+
| gaxpy | X | |
+----------------------+-------------+---------+
| get_scalar | X | |
+----------------------+-------------+---------+
| has_nan | X | |
+----------------------+-------------+---------+
| init (splitlink) | X | |
+----------------------+-------------+---------+
| is_eye | X | |
+----------------------+-------------+---------+
| is_set | X | |
+----------------------+-------------+---------+
| kron | X | |
+----------------------+-------------+---------+
| maxlineardim | X | |
+----------------------+-------------+---------+
| maxvalue | X | |
+----------------------+-------------+---------+
| norm | X | |
+----------------------+-------------+---------+
| par_recv | | X |
+----------------------+-------------+---------+
| par_send | | X |
+----------------------+-------------+---------+
| permute | X | |
+----------------------+-------------+---------+
| perturb | X | |
+----------------------+-------------+---------+
| pointto | X | |
+----------------------+-------------+---------+
| print | X | |
+----------------------+-------------+---------+
| project | X | |
+----------------------+-------------+---------+
| qmattomat | X | |
+----------------------+-------------+---------+
| randomize | X | |
+----------------------+-------------+---------+
| rank | X | |
+----------------------+-------------+---------+
| read | X | |
+----------------------+-------------+---------+
| scale | X | |
+----------------------+-------------+---------+
| set_hash | X | |
+----------------------+-------------+---------+
| split | X | |
+----------------------+-------------+---------+
| trace | X | |
+----------------------+-------------+---------+
| transposed | X | |
+----------------------+-------------+---------+
| write | X | |
+----------------------+-------------+---------+
| write_as_matrix | X | |
+----------------------+-------------+---------+
"""
[docs]def addProject_tensor():
"""
fortran-subroutine - March 2016 (updated dj)
Apply the projector P to A and add to B, where P is defined as
:math:`P = 1-\sum_{\\alpha}|psiProjs_{\\alpha}><psiProjs_{\\alpha}|`.
**Arguments**
Bb : TYPE(tensor), inout
Add projection of A to this tensor.
Aa : TYPE(tensor), in
Get projection of this tensor and add to Bb
PsiProjs : TYPE(tensor)(\*), in
Array of tensors defining the projector.
**Details**
Used in orthogonalizing an MPS against another set of MPSs.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine addProject_tensor(Bb, Aa, PsiProjs, errst)
type(tensor), intent(inout) :: Bb
type(tensor), intent(in) :: Aa
type(tensor), pointer, intent(in) :: PsiProjs(:)
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: jj
! overlap of tensors (dot product)
real(KIND=rKind) :: proj
!if(present(errst)) errst = 0
! For done we can use underlying BLAS subroutines
call gaxpy(Bb, done, Aa)
do jj = 1, size(PsiProjs)
proj = - dot(PsiProjs(jj), Aa)
call gaxpy(Bb, proj, PsiProjs(jj))
end do
end subroutine addProject_tensor
"""
return
[docs]def addProject_tensorc():
"""
fortran-subroutine - March 2016 (updated dj)
Apply the projector P to A and add to B, where P is defined as
:math:`P = 1-\sum_{\\alpha}|psiProjs_{\\alpha}><psiProjs_{\\alpha}|`.
**Arguments**
Bb : TYPE(tensorc), inout
Add projection of A to this tensor.
Aa : TYPE(tensorc), in
Get projection of this tensor and add to Bb
PsiProjs : TYPE(tensorc)(\*), in
Array of tensors defining the projector.
**Details**
Used in orthogonalizing an MPS against another set of MPSs.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine addProject_tensorc(Bb, Aa, PsiProjs, errst)
type(tensorc), intent(inout) :: Bb
type(tensorc), intent(in) :: Aa
type(tensorc), pointer, intent(in) :: PsiProjs(:)
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: jj
! overlap of tensors (dot product)
complex(KIND=rKind) :: proj
!if(present(errst)) errst = 0
! For zone we can use underlying BLAS subroutines
call gaxpy(Bb, zone, Aa)
do jj = 1, size(PsiProjs)
proj = - dot(PsiProjs(jj), Aa)
call gaxpy(Bb, proj, PsiProjs(jj))
end do
end subroutine addProject_tensorc
"""
return
[docs]def add_dummylink_tensor():
"""
fortran-subroutine - October 2017 (dj)
Add a dummy link to a tensor.
**Arguments**
Tens : TYPE(tensor), inout
Add a dummy link of dimenions one to the tensor.
idx : INTEGER, in
Position of the dummy link in the modified tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine add_dummylink_tensor(Tens, idx, errst)
type(tensor), intent(inout) :: Tens
integer, intent(in) :: idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! new dimension of the tensor
integer, dimension(:), allocatable :: newdl
allocate(newdl(Tens%rank + 1))
! Increase rank
Tens%rank = Tens%rank + 1
! Modify dl
newdl(:idx - 1) = Tens%dl(:idx - 1)
newdl(idx) = 1
newdl(idx + 1:) = Tens%dl(idx:)
deallocate(Tens%dl)
allocate(Tens%dl(Tens%rank))
Tens%dl = newdl
! Modify mdl
newdl(:idx - 1) = Tens%mdl(:idx - 1)
newdl(idx) = 1
newdl(idx + 1:) = Tens%mdl(idx:)
deallocate(Tens%mdl)
allocate(Tens%mdl(Tens%rank))
Tens%mdl = newdl
! Modify Tens%idx
newdl(:idx - 1) = Tens%idx(:idx - 1)
newdl(idx) = idx
newdl(idx + 1:) = Tens%idx(idx:) + 1
deallocate(Tens%idx)
allocate(Tens%idx(Tens%rank))
Tens%idx = newdl
deallocate(newdl)
end subroutine add_dummylink_tensor
"""
return
[docs]def add_dummylink_tensorc():
"""
fortran-subroutine - October 2017 (dj)
Add a dummy link to a tensor.
**Arguments**
Tens : TYPE(tensorc), inout
Add a dummy link of dimenions one to the tensor.
idx : INTEGER, in
Position of the dummy link in the modified tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine add_dummylink_tensorc(Tens, idx, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! new dimension of the tensor
integer, dimension(:), allocatable :: newdl
allocate(newdl(Tens%rank + 1))
! Increase rank
Tens%rank = Tens%rank + 1
! Modify dl
newdl(:idx - 1) = Tens%dl(:idx - 1)
newdl(idx) = 1
newdl(idx + 1:) = Tens%dl(idx:)
deallocate(Tens%dl)
allocate(Tens%dl(Tens%rank))
Tens%dl = newdl
! Modify mdl
newdl(:idx - 1) = Tens%mdl(:idx - 1)
newdl(idx) = 1
newdl(idx + 1:) = Tens%mdl(idx:)
deallocate(Tens%mdl)
allocate(Tens%mdl(Tens%rank))
Tens%mdl = newdl
! Modify Tens%idx
newdl(:idx - 1) = Tens%idx(:idx - 1)
newdl(idx) = idx
newdl(idx + 1:) = Tens%idx(idx:) + 1
deallocate(Tens%idx)
allocate(Tens%idx(Tens%rank))
Tens%idx = newdl
deallocate(newdl)
end subroutine add_dummylink_tensorc
"""
return
[docs]def check_qnum_tensdim_tensor():
"""
fortran-subroutine - August 2018 (dj)
Dummy routine to comply with qTensors.
**Arguments**
Tens : TYPE(Qtensor), in
Tensor to be checked. Nothing to check.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine check_qnum_tensdim_tensor(Tens, errst)
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii, jj
!if(present(errst)) errst = 0
end subroutine check_qnum_tensdim_tensor
"""
return
[docs]def check_qnum_tensdim_tensorc():
"""
fortran-subroutine - August 2018 (dj)
Dummy routine to comply with qTensors.
**Arguments**
Tens : TYPE(Qtensorc), in
Tensor to be checked. Nothing to check.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine check_qnum_tensdim_tensorc(Tens, errst)
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii, jj
!if(present(errst)) errst = 0
end subroutine check_qnum_tensdim_tensorc
"""
return
[docs]def conj_tensor():
"""
fortran-subroutine - October 2016 (dj)
Complex-conjugated elements in the tensor.
**Arguments**
Tens : TYPE(tensor), inout
Take the complex conjugate of each element
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine conj_tensor(Tens, errst)
type(tensor), intent(inout) :: Tens
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!Tens%elem(:Tens%dim) = conjg(Tens%elem(:Tens%dim))
end subroutine conj_tensor
"""
return
[docs]def conj_tensorc():
"""
fortran-subroutine - October 2016 (dj)
Complex-conjugated elements in the tensor.
**Arguments**
Tens : TYPE(tensorc), inout
Take the complex conjugate of each element
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine conj_tensorc(Tens, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
Tens%elem(:Tens%dim) = conjg(Tens%elem(:Tens%dim))
end subroutine conj_tensorc
"""
return
[docs]def copy_tensor_tensor():
"""
fortran-subroutine - October 2016 (dj)
Copy a tensor.
**Arguments**
Tensout : TYPE(tensor), out
Store copy of Tensin in this tensor.
Tensin : TYPE(tensor), in
Copy this tensor to a new tensor Tensout.
scalar : real, OPTIONAL, in
Multiply input tensor with scalar during copying replacing
a copy + scale action.
trans : CHARACTER, OPTIONAL, in
The following transformation can be applied: complex conjugate ('C'),
transposition ('T', simple transposition, no permutation), conjugate
transposed ('H', simple transposition). No transformation is 'N'.
Default to 'N'.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensor_tensor(Tensout, Tensin, scalar, trans, &
errst)
type(tensor), intent(out) :: Tensout
type(tensor), intent(in) :: Tensin
real(KIND=rKind), intent(in), optional :: scalar
character, intent(in), optional :: trans
integer, intent(out), optional :: errst
!if(Tensin%rank == -1) then
! errst = raise_error('copy_tensor_tensor: '//&
! 'Tensin empty.', 99, 'Tensors_include.f90:263', &
! errst=errst)
! return
!end if
call create(Tensout, Tensin%dl, errst=errst)
!if(prop_error('copy_tensor_tensor: create failed.', &
! 'Tensors_include.f90:270', errst=errst)) return
if(present(scalar)) then
Tensout%elem(:Tensin%dim) = scalar * Tensin%elem(:Tensin%dim)
else
Tensout%elem(:Tensin%dim) = Tensin%elem(:Tensin%dim)
end if
! No transformation
if(.not. present(trans)) return
! Check which transformation is applied
if(trans == 'C') then
call conj(Tensout, errst)
!if(prop_error('copy_tensor_tensor:'//&
! ' conjg failed.', 'Tensors_include.f90:285', &
! errst=errst)) return
elseif(trans == 'T') then
call transposed(Tensout, errst=errst)
!if(prop_error('copy_tensor_tensor:'//&
! ' transpose failed.', 'Tensors_include.f90:291', &
! errst=errst)) return
elseif(trans == 'H') then
call dagger(Tensout, errst=errst)
!if(prop_error('copy_tensor_tensor:'//&
! 'dagger failed.', 'Tensors_include.f90:297', &
! errst=errst)) return
elseif(trans /= 'N') then
!errst = raise_error('copy_tensor_tensor:'//&
! ' unknown argument trans.', 4, &
! 'Tensors_include.f90:303', errst=errst)
end if
end subroutine copy_tensor_tensor
"""
return
[docs]def copy_tensorc_tensor():
"""
fortran-subroutine - October 2016 (dj)
Copy a tensor.
**Arguments**
Tensout : TYPE(tensorc), out
Store copy of Tensin in this tensor.
Tensin : TYPE(tensor), in
Copy this tensor to a new tensor Tensout.
scalar : complex, OPTIONAL, in
Multiply input tensor with scalar during copying replacing
a copy + scale action.
trans : CHARACTER, OPTIONAL, in
The following transformation can be applied: complex conjugate ('C'),
transposition ('T', simple transposition, no permutation), conjugate
transposed ('H', simple transposition). No transformation is 'N'.
Default to 'N'.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensorc_tensor(Tensout, Tensin, scalar, trans, &
errst)
type(tensorc), intent(out) :: Tensout
type(tensor), intent(in) :: Tensin
complex(KIND=rKind), intent(in), optional :: scalar
character, intent(in), optional :: trans
integer, intent(out), optional :: errst
!if(Tensin%rank == -1) then
! errst = raise_error('copy_tensorc_tensor: '//&
! 'Tensin empty.', 99, 'Tensors_include.f90:263', &
! errst=errst)
! return
!end if
call create(Tensout, Tensin%dl, errst=errst)
!if(prop_error('copy_tensorc_tensor: create failed.', &
! 'Tensors_include.f90:270', errst=errst)) return
if(present(scalar)) then
Tensout%elem(:Tensin%dim) = scalar * Tensin%elem(:Tensin%dim)
else
Tensout%elem(:Tensin%dim) = Tensin%elem(:Tensin%dim)
end if
! No transformation
if(.not. present(trans)) return
! Check which transformation is applied
if(trans == 'C') then
call conj(Tensout, errst)
!if(prop_error('copy_tensorc_tensor:'//&
! ' conjg failed.', 'Tensors_include.f90:285', &
! errst=errst)) return
elseif(trans == 'T') then
call transposed(Tensout, errst=errst)
!if(prop_error('copy_tensorc_tensor:'//&
! ' transpose failed.', 'Tensors_include.f90:291', &
! errst=errst)) return
elseif(trans == 'H') then
call dagger(Tensout, errst=errst)
!if(prop_error('copy_tensorc_tensor:'//&
! 'dagger failed.', 'Tensors_include.f90:297', &
! errst=errst)) return
elseif(trans /= 'N') then
!errst = raise_error('copy_tensorc_tensor:'//&
! ' unknown argument trans.', 4, &
! 'Tensors_include.f90:303', errst=errst)
end if
end subroutine copy_tensorc_tensor
"""
return
[docs]def copy_tensorc_tensorc():
"""
fortran-subroutine - October 2016 (dj)
Copy a tensor.
**Arguments**
Tensout : TYPE(tensorc), out
Store copy of Tensin in this tensor.
Tensin : TYPE(tensorc), in
Copy this tensor to a new tensor Tensout.
scalar : complex, OPTIONAL, in
Multiply input tensor with scalar during copying replacing
a copy + scale action.
trans : CHARACTER, OPTIONAL, in
The following transformation can be applied: complex conjugate ('C'),
transposition ('T', simple transposition, no permutation), conjugate
transposed ('H', simple transposition). No transformation is 'N'.
Default to 'N'.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensorc_tensorc(Tensout, Tensin, scalar, trans, &
errst)
type(tensorc), intent(out) :: Tensout
type(tensorc), intent(in) :: Tensin
complex(KIND=rKind), intent(in), optional :: scalar
character, intent(in), optional :: trans
integer, intent(out), optional :: errst
!if(Tensin%rank == -1) then
! errst = raise_error('copy_tensorc_tensorc: '//&
! 'Tensin empty.', 99, 'Tensors_include.f90:263', &
! errst=errst)
! return
!end if
call create(Tensout, Tensin%dl, errst=errst)
!if(prop_error('copy_tensorc_tensorc: create failed.', &
! 'Tensors_include.f90:270', errst=errst)) return
if(present(scalar)) then
Tensout%elem(:Tensin%dim) = scalar * Tensin%elem(:Tensin%dim)
else
Tensout%elem(:Tensin%dim) = Tensin%elem(:Tensin%dim)
end if
! No transformation
if(.not. present(trans)) return
! Check which transformation is applied
if(trans == 'C') then
call conj(Tensout, errst)
!if(prop_error('copy_tensorc_tensorc:'//&
! ' conjg failed.', 'Tensors_include.f90:285', &
! errst=errst)) return
elseif(trans == 'T') then
call transposed(Tensout, errst=errst)
!if(prop_error('copy_tensorc_tensorc:'//&
! ' transpose failed.', 'Tensors_include.f90:291', &
! errst=errst)) return
elseif(trans == 'H') then
call dagger(Tensout, errst=errst)
!if(prop_error('copy_tensorc_tensorc:'//&
! 'dagger failed.', 'Tensors_include.f90:297', &
! errst=errst)) return
elseif(trans /= 'N') then
!errst = raise_error('copy_tensorc_tensorc:'//&
! ' unknown argument trans.', 4, &
! 'Tensors_include.f90:303', errst=errst)
end if
end subroutine copy_tensorc_tensorc
"""
return
[docs]def copy_tensorlist():
"""
fortran-subroutine - ?? ()
Copy a tensor list.
**Arguments**
Objout : TYPE(tensorlist), out
Copy the entries of `Objin` here.
Objin : TYPE(tensorlist), in
source for copy into `Dest`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensorlist(Objout, Objin)
type(tensorlist), intent(out) :: Objout
type(tensorlist), intent(in) :: Objin
! Local variables
! ---------------
! for looping
integer :: ii
allocate(Objout%Li(size(Objin%Li)))
do ii = 1, size(Objin%Li)
call copy(Objout%Li(ii), Objin%Li(ii))
end do
end subroutine copy_tensorlist
"""
return
[docs]def copy_tensorlistc():
"""
fortran-subroutine - ?? ()
Copy a tensorc list.
**Arguments**
Objout : TYPE(tensorlistc), out
Copy the entries of `Objin` here.
Objin : TYPE(tensorlistc), in
source for copy into `Dest`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensorlistc(Objout, Objin)
type(tensorlistc), intent(out) :: Objout
type(tensorlistc), intent(in) :: Objin
! Local variables
! ---------------
! for looping
integer :: ii
allocate(Objout%Li(size(Objin%Li)))
do ii = 1, size(Objin%Li)
call copy(Objout%Li(ii), Objin%Li(ii))
end do
end subroutine copy_tensorlistc
"""
return
[docs]def copy_tensorlistarray():
"""
fortran-subroutine - ?? ()
Copy an array of tensor lists.
**Arguments**
Objout : TYPE(tensorlist)(\*), out
Copy the entries of `Objin` here.
Objin : TYPE(tensorlist)(\*), in
source for copy into `Objout`
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensorlistarray(Objout, Objin)
type(tensorlist), pointer, intent(out) :: Objout(:)
type(tensorlist), pointer, intent(in) :: Objin(:)
! Local variables
! ---------------
! for looping
integer :: ii
allocate(Objout(size(Objin)))
do ii = 1, size(Objin)
call copy(Objout(ii), Objin(ii))
end do
end subroutine copy_tensorlistarray
"""
return
[docs]def copy_tensorlistarrayc():
"""
fortran-subroutine - ?? ()
Copy an array of tensorc lists.
**Arguments**
Objout : TYPE(tensorlistc)(\*), out
Copy the entries of `Objin` here.
Objin : TYPE(tensorlistc)(\*), in
source for copy into `Objout`
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_tensorlistarrayc(Objout, Objin)
type(tensorlistc), pointer, intent(out) :: Objout(:)
type(tensorlistc), pointer, intent(in) :: Objin(:)
! Local variables
! ---------------
! for looping
integer :: ii
allocate(Objout(size(Objin)))
do ii = 1, size(Objin)
call copy(Objout(ii), Objin(ii))
end do
end subroutine copy_tensorlistarrayc
"""
return
[docs]def create_tensor():
"""
fortran-subroutine - July 2016 (dj)
Create a tensor of rank-k.
**Arguments**
Tens : TYPE(tensor), out
Allocate an rank-k tensor with certain dimensions.
dl : INTEGER, in
Dimension of the vector to be allocated.
init : CHARACTER, OPTIONAL, in
Options for initialization of the array. These are initialization
with zeros ('0'), ones ('1'), no initialization ('N') and random
numbers ('R').
Default to 'N' (if not present).
mdl : INTEGER(\*), OPTIONAL, in
Maximal dimension of the tensor during its whole use in the
algorithm. Allows to avoid deallocation/allocation steps.
Default to dl.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_tensor(Tens, dl, init, mdl, errst)
type(tensor), intent(inout) :: Tens
integer, dimension(:), intent(in) :: dl
character, intent(in), optional :: init
integer, dimension(:), intent(in), optional :: mdl
integer, intent(out), optional :: errst
! Local variable
! --------------
! for looping
integer :: ii
!if(present(errst)) errst = 0
if(associated(Tens%elem)) then
call destroy(Tens, errst=errst)
!if(prop_error('create_tensor: dealloc failed.', &
! 'Tensors_include.f90:459', errst=errst)) return
end if
Tens%rank = size(dl)
if(Tens%rank == 0) then
! Scalar
Tens%dim = 1
allocate(Tens%elem(1))
else
! Tensor
allocate(Tens%dl(Tens%rank))
Tens%dl = dl
Tens%dim = product(dl)
allocate(Tens%mdl(Tens%rank))
if(present(mdl)) then
Tens%mdl = mdl
else
Tens%mdl = dl
end if
allocate(Tens%idx(Tens%rank))
Tens%idx = [(ii, ii = 1, Tens%rank)]
allocate(Tens%elem(product(Tens%mdl)))
end if
! No initialization
if(.not. present(init)) then
return
end if
! Initialize with given option
if(init == '0') then
Tens%elem(:Tens%dim) = dzero
elseif(init == '1') then
Tens%elem(:Tens%dim) = done
elseif(init == 'R') then
call randomize(Tens, errst=errst)
!if(prop_error('create_tensor: randomize failed.', &
! 'Tensors_include.f90:499', errst=errst)) return
elseif(init /= 'N') then
!errst = raise_error('create_tensor: unbknown option', &
! 4, 'Tensors_include.f90:502', errst=errst)
end if
end subroutine create_tensor
"""
return
[docs]def create_tensorc():
"""
fortran-subroutine - July 2016 (dj)
Create a tensor of rank-k.
**Arguments**
Tens : TYPE(tensorc), out
Allocate an rank-k tensor with certain dimensions.
dl : INTEGER, in
Dimension of the vector to be allocated.
init : CHARACTER, OPTIONAL, in
Options for initialization of the array. These are initialization
with zeros ('0'), ones ('1'), no initialization ('N') and random
numbers ('R').
Default to 'N' (if not present).
mdl : INTEGER(\*), OPTIONAL, in
Maximal dimension of the tensor during its whole use in the
algorithm. Allows to avoid deallocation/allocation steps.
Default to dl.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_tensorc(Tens, dl, init, mdl, errst)
type(tensorc), intent(inout) :: Tens
integer, dimension(:), intent(in) :: dl
character, intent(in), optional :: init
integer, dimension(:), intent(in), optional :: mdl
integer, intent(out), optional :: errst
! Local variable
! --------------
! for looping
integer :: ii
!if(present(errst)) errst = 0
if(associated(Tens%elem)) then
call destroy(Tens, errst=errst)
!if(prop_error('create_tensorc: dealloc failed.', &
! 'Tensors_include.f90:459', errst=errst)) return
end if
Tens%rank = size(dl)
if(Tens%rank == 0) then
! Scalar
Tens%dim = 1
allocate(Tens%elem(1))
else
! Tensor
allocate(Tens%dl(Tens%rank))
Tens%dl = dl
Tens%dim = product(dl)
allocate(Tens%mdl(Tens%rank))
if(present(mdl)) then
Tens%mdl = mdl
else
Tens%mdl = dl
end if
allocate(Tens%idx(Tens%rank))
Tens%idx = [(ii, ii = 1, Tens%rank)]
allocate(Tens%elem(product(Tens%mdl)))
end if
! No initialization
if(.not. present(init)) then
return
end if
! Initialize with given option
if(init == '0') then
Tens%elem(:Tens%dim) = zzero
elseif(init == '1') then
Tens%elem(:Tens%dim) = zone
elseif(init == 'R') then
call randomize(Tens, errst=errst)
!if(prop_error('create_tensorc: randomize failed.', &
! 'Tensors_include.f90:499', errst=errst)) return
elseif(init /= 'N') then
!errst = raise_error('create_tensorc: unbknown option', &
! 4, 'Tensors_include.f90:502', errst=errst)
end if
end subroutine create_tensorc
"""
return
[docs]def create_id_tensor():
"""
fortran-subroutine - October 2016 (updated dj)
Create an identity matrix (or delta function, rank-2 tensor).
**Arguments**
Idmat : TYPE(tensor), out
Create and identity matrix of dimension `dim x dim`.
Tens : TYPE(tensor), in
Needed to extract dimension for identity matrix. Has to be rank-3.
side : INTEGER, in
Either `1` or `3` depending if identity matrix should be created
as LeftDelta (`1`, first index) or RightDelta (`3`, last/third
index).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_id_tensor(Idmat, Tens, side, errst)
type(tensor), intent(out) :: Idmat
type(tensor), intent(in) :: Tens
integer, intent(in) :: side
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! storing dimension
integer :: dim
!if(present(errst)) errst = 0
!if(Tens%rank /= 3) then
! errst = raise_error(&
! 'create_id_tensor: need rank-3 tensor.', &
! 2, 'Tensors_include.f90:561', errst=errst)
! return
!end if
!if(.not. ((side == 1) .or. (side == 3))) then
! errst = raise_error(&
! 'create_id_tensor: side must be 1 or 3.', &
! 4, 'Tensors_include.f90:568', errst=errst)
! return
!end if
dim = Tens%dl(side)
call create(Idmat, [dim, dim], init='0', errst=errst)
!if(prop_error('create_id_tensor: create failed.', &
! 'Tensors_include.f90:575', errst=errst)) return
do ii = 1, dim
Idmat%elem(ii + (ii - 1) * dim) = done
end do
end subroutine create_id_tensor
"""
return
[docs]def create_id_tensorc():
"""
fortran-subroutine - October 2016 (updated dj)
Create an identity matrix (or delta function, rank-2 tensor).
**Arguments**
Idmat : TYPE(tensorc), out
Create and identity matrix of dimension `dim x dim`.
Tens : TYPE(tensorc), in
Needed to extract dimension for identity matrix. Has to be rank-3.
side : INTEGER, in
Either `1` or `3` depending if identity matrix should be created
as LeftDelta (`1`, first index) or RightDelta (`3`, last/third
index).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_id_tensorc(Idmat, Tens, side, errst)
type(tensorc), intent(out) :: Idmat
type(tensorc), intent(in) :: Tens
integer, intent(in) :: side
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! storing dimension
integer :: dim
!if(present(errst)) errst = 0
!if(Tens%rank /= 3) then
! errst = raise_error(&
! 'create_id_tensorc: need rank-3 tensor.', &
! 2, 'Tensors_include.f90:561', errst=errst)
! return
!end if
!if(.not. ((side == 1) .or. (side == 3))) then
! errst = raise_error(&
! 'create_id_tensorc: side must be 1 or 3.', &
! 4, 'Tensors_include.f90:568', errst=errst)
! return
!end if
dim = Tens%dl(side)
call create(Idmat, [dim, dim], init='0', errst=errst)
!if(prop_error('create_id_tensorc: create failed.', &
! 'Tensors_include.f90:575', errst=errst)) return
do ii = 1, dim
Idmat%elem(ii + (ii - 1) * dim) = zone
end do
end subroutine create_id_tensorc
"""
return
[docs]def create_diag_tensor_tensor():
"""
fortran-subroutine - June 2018 (dj)
Create a rank-2 tensor with the given values on the diagonal.
**Arguments**
Dmat : TYPE(tensor), inout
On exit, diagonal rank-2 tensor.
Dvals : TYPE(tensor), in
Diagonal entries for the rank-2 tensor as rank-1 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_diag_tensor_tensor(Dmat, Dvals, errst)
type(tensor), intent(inout) :: Dmat
type(tensor), intent(in) :: Dvals
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! index for next diagonal element
integer :: jj
! storing dimension
integer :: dim
!if(present(errst)) errst = 0
dim = Dvals%dl(1)
call create(Dmat, [dim, dim], init='0', errst=errst)
!if(prop_error('create_diag_tensor: create failed.', &
! 'Tensors_include.f90:631', errst=errst)) return
jj = 1
do ii = 1, dim
Dmat%elem(jj) = Dvals%elem(ii)
jj = jj + dim + 1
end do
end subroutine create_diag_tensor_tensor
"""
return
[docs]def create_diag_tensor_tensorc():
"""
fortran-subroutine - June 2018 (dj)
Create a rank-2 tensor with the given values on the diagonal.
**Arguments**
Dmat : TYPE(tensorc), inout
On exit, diagonal rank-2 tensor.
Dvals : TYPE(tensor), in
Diagonal entries for the rank-2 tensor as rank-1 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_diag_tensor_tensorc(Dmat, Dvals, errst)
type(tensorc), intent(inout) :: Dmat
type(tensor), intent(in) :: Dvals
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! index for next diagonal element
integer :: jj
! storing dimension
integer :: dim
!if(present(errst)) errst = 0
dim = Dvals%dl(1)
call create(Dmat, [dim, dim], init='0', errst=errst)
!if(prop_error('create_diag_tensorc: create failed.', &
! 'Tensors_include.f90:631', errst=errst)) return
jj = 1
do ii = 1, dim
Dmat%elem(jj) = Dvals%elem(ii)
jj = jj + dim + 1
end do
end subroutine create_diag_tensor_tensorc
"""
return
[docs]def create_diag_tensorc_tensorc():
"""
fortran-subroutine - June 2018 (dj)
Create a rank-2 tensor with the given values on the diagonal.
**Arguments**
Dmat : TYPE(tensorc), inout
On exit, diagonal rank-2 tensor.
Dvals : TYPE(tensorc), in
Diagonal entries for the rank-2 tensor as rank-1 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_diag_tensorc_tensorc(Dmat, Dvals, errst)
type(tensorc), intent(inout) :: Dmat
type(tensorc), intent(in) :: Dvals
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! index for next diagonal element
integer :: jj
! storing dimension
integer :: dim
!if(present(errst)) errst = 0
dim = Dvals%dl(1)
call create(Dmat, [dim, dim], init='0', errst=errst)
!if(prop_error('create_diag_tensorc: create failed.', &
! 'Tensors_include.f90:631', errst=errst)) return
jj = 1
do ii = 1, dim
Dmat%elem(jj) = Dvals%elem(ii)
jj = jj + dim + 1
end do
end subroutine create_diag_tensorc_tensorc
"""
return
[docs]def create_tensorlist():
"""
fortran-subroutine - January 2016 (update dj)
Create a tensor list of length d1.
**Arguments**
Mats : TYPE(tensorlist), inout
Allocate an array of size d1.
d1 : INTEGER, in
Length of the matrix list to be allocated.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_tensorlist(Mats, d1)
type(tensorlist), intent(inout) :: Mats
integer, intent(in) :: d1
allocate(Mats%Li(d1))
end subroutine create_tensorlist
"""
return
[docs]def create_tensorlistc():
"""
fortran-subroutine - January 2016 (update dj)
Create a tensorc list of length d1.
**Arguments**
Mats : TYPE(tensorlistc), inout
Allocate an array of size d1.
d1 : INTEGER, in
Length of the matrix list to be allocated.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_tensorlistc(Mats, d1)
type(tensorlistc), intent(inout) :: Mats
integer, intent(in) :: d1
allocate(Mats%Li(d1))
end subroutine create_tensorlistc
"""
return
[docs]def dagger_tensor():
"""
fortran-subroutine - October 2016 (dj)
Taking complex conjugate followed by transposition of indices or
permutation.
**Arguments**
Tens : TYPE(tensor), inout
Save a transposition/permutation on the indices of this tensor.
Complex conjugated values are taken during this subroutine.
perm : INTEGER(\*), OPTIONAL, in
permutation array has length equal to the rank of the tensor with
unique entries 1 to rank.
Default to rank, rank - 1, ..., 2, 1 (transpose)
doperm : LOGICAL, OPTIONAL, in
If true, then permute immediately. If false or not present,
permutation is only indicated in the indices.
**Details**
For details of permutation look into
Tensors.f90:transpose_tensor
(template defined in Tensors_include.f90)
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine dagger_tensor(Tens, perm, doperm, errst)
type(tensor), intent(inout) :: Tens
integer, dimension(:), intent(in), optional :: perm
logical, intent(in), optional :: doperm
integer, intent(out), optional :: errst
call conj(Tens, errst=errst)
!if(prop_error('dagger_tensor: conjg failed.', &
! 'Tensors_include.f90:726', errst=errst)) return
call transposed(Tens, perm, doperm, errst=errst)
!if(prop_error('dagger_tensor: transpose failed.', &
! 'Tensors_include.f90:730', errst=errst)) return
end subroutine dagger_tensor
"""
return
[docs]def destroy_tensor():
"""
fortran-subroutine - July 2016 (dj)
Destroy a tensor.
**Arguments**
Tens : TYPE(tensor), inout
Deallocate the array representing the tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine destroy_tensor(Tens, errst)
type(tensor), intent(inout) :: Tens
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(.not. associated(Tens%elem)) then
! errst = raise_error('destroy_tensor: elem not '// &
! 'allocated.', 1, 'Tensors_include.f90:758', errst=errst)
! return
!end if
if(Tens%rank == 0) then
! Scalar
deallocate(Tens%elem)
Tens%rank = -1
Tens%dim = 0
else
! Tensor
deallocate(Tens%dl, Tens%mdl, Tens%idx, Tens%elem)
Tens%rank = -1
Tens%dim = 0
end if
end subroutine destroy_tensor
"""
return
[docs]def dagger_tensorc():
"""
fortran-subroutine - October 2016 (dj)
Taking complex conjugate followed by transposition of indices or
permutation.
**Arguments**
Tens : TYPE(tensorc), inout
Save a transposition/permutation on the indices of this tensor.
Complex conjugated values are taken during this subroutine.
perm : INTEGER(\*), OPTIONAL, in
permutation array has length equal to the rank of the tensor with
unique entries 1 to rank.
Default to rank, rank - 1, ..., 2, 1 (transpose)
doperm : LOGICAL, OPTIONAL, in
If true, then permute immediately. If false or not present,
permutation is only indicated in the indices.
**Details**
For details of permutation look into
Tensors.f90:transpose_tensorc
(template defined in Tensors_include.f90)
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine dagger_tensorc(Tens, perm, doperm, errst)
type(tensorc), intent(inout) :: Tens
integer, dimension(:), intent(in), optional :: perm
logical, intent(in), optional :: doperm
integer, intent(out), optional :: errst
call conj(Tens, errst=errst)
!if(prop_error('dagger_tensorc: conjg failed.', &
! 'Tensors_include.f90:726', errst=errst)) return
call transposed(Tens, perm, doperm, errst=errst)
!if(prop_error('dagger_tensorc: transpose failed.', &
! 'Tensors_include.f90:730', errst=errst)) return
end subroutine dagger_tensorc
"""
return
[docs]def destroy_tensorc():
"""
fortran-subroutine - July 2016 (dj)
Destroy a tensor.
**Arguments**
Tens : TYPE(tensorc), inout
Deallocate the array representing the tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine destroy_tensorc(Tens, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(.not. associated(Tens%elem)) then
! errst = raise_error('destroy_tensorc: elem not '// &
! 'allocated.', 1, 'Tensors_include.f90:758', errst=errst)
! return
!end if
if(Tens%rank == 0) then
! Scalar
deallocate(Tens%elem)
Tens%rank = -1
Tens%dim = 0
else
! Tensor
deallocate(Tens%dl, Tens%mdl, Tens%idx, Tens%elem)
Tens%rank = -1
Tens%dim = 0
end if
end subroutine destroy_tensorc
"""
return
[docs]def destroy_tensorlist():
"""
fortran-subroutine - January 2016 (update dj)
Destroy a list of tensor.
**Arguments**
Mats : TYPE(tensorlist), inout
Deallocate the array of matrixs (including each matrix).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine destroy_tensorlist(Mats)
type(tensorlist), intent(inout) :: Mats
! Local variables
! ---------------
! for looping
integer :: ii
do ii = 1, size(Mats%Li)
call destroy(Mats%Li(ii))
end do
deallocate(Mats%Li)
end subroutine destroy_tensorlist
"""
return
[docs]def destroy_tensorlistc():
"""
fortran-subroutine - January 2016 (update dj)
Destroy a list of tensorc.
**Arguments**
Mats : TYPE(tensorlistc), inout
Deallocate the array of matrixs (including each matrix).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine destroy_tensorlistc(Mats)
type(tensorlistc), intent(inout) :: Mats
! Local variables
! ---------------
! for looping
integer :: ii
do ii = 1, size(Mats%Li)
call destroy(Mats%Li(ii))
end do
deallocate(Mats%Li)
end subroutine destroy_tensorlistc
"""
return
[docs]def destroy_imap():
"""
fortran-subroutine - September 2017 (dj)
Destroy the imap if allocated.
**Arguments**
Imapper : TYPE(imap), inout
Deallocated on exit. It checks if allocated first.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine destroy_imap(Imapper)
type(imap), intent(inout) :: Imapper
! No local variables
! ------------------
if(allocated(Imapper%start)) then
deallocate(Imapper%start, Imapper%dd, Imapper%hashes)
end if
end subroutine destroy_imap
"""
return
[docs]def dot_tensor():
"""
fortran-function - October 2016 (dj)
Compute dot product <A, B> for tensors.
**Arguments**
Tensa : TYPE(tensor), out
First tensor in the dot product (taken with complex conjugated values
for complex arrays).
Tensb : TYPE(tensor), in
Second tensor in the dot product. Has to be the same size as first
tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function dot_tensor(Tensa, Tensb, errst) result(dot)
type(tensor), intent(in) :: Tensa, Tensb
integer, intent(out), optional :: errst
real(KIND=rKind) :: dot
! External functions
! ------------------
! goes into DDOT,ZDOTC
real(KIND=rKind), external :: DDOT
!if(present(errst)) errst = 0
!if(Tensa%rank /= Tensb%rank) then
! errst = raise_error('dot_tensor: rank mismatch.', &
! 2, 'Tensors_include.f90:888', errst=errst)
! return
!end if
!if((Tensa%rank /= 0) .and. any(Tensa%dl /= Tensb%dl)) then
! errst = raise_error('dot_tensor: dimension mismatch.', &
! 3, 'Tensors_include.f90:894', errst=errst)
! return
!end if
dot = DDOT(Tensa%dim, Tensa%elem, 1, Tensb%elem, 1)
!dot = macdot(Tensa%dim, Tensa%elem, Tensb%elem)
end function dot_tensor
"""
return
[docs]def dot_tensorc():
"""
fortran-function - October 2016 (dj)
Compute dot product <A, B> for tensors.
**Arguments**
Tensa : TYPE(tensorc), out
First tensor in the dot product (taken with complex conjugated values
for complex arrays).
Tensb : TYPE(tensorc), in
Second tensor in the dot product. Has to be the same size as first
tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function dot_tensorc(Tensa, Tensb, errst) result(dot)
type(tensorc), intent(in) :: Tensa, Tensb
integer, intent(out), optional :: errst
complex(KIND=rKind) :: dot
! External functions
! ------------------
! goes into DDOT,ZDOTC
complex(KIND=rKind), external :: ZDOTC
!if(present(errst)) errst = 0
!if(Tensa%rank /= Tensb%rank) then
! errst = raise_error('dot_tensorc: rank mismatch.', &
! 2, 'Tensors_include.f90:888', errst=errst)
! return
!end if
!if((Tensa%rank /= 0) .and. any(Tensa%dl /= Tensb%dl)) then
! errst = raise_error('dot_tensorc: dimension mismatch.', &
! 3, 'Tensors_include.f90:894', errst=errst)
! return
!end if
dot = ZDOTC(Tensa%dim, Tensa%elem, 1, Tensb%elem, 1)
!dot = macdot(Tensa%dim, Tensa%elem, Tensb%elem)
end function dot_tensorc
"""
return
[docs]def macdot_real():
"""
fortran-function - July 2018 (dj)
Dot product implemented to avoid faulty ZDOTC implementation of Mac OS X
accelerate package.
**Arguments**
dim : INTEGER, in
Dimension of the dot product.
elema : real(\*), in
First argument in dot product. Complex-conjugate taken
for complex-valued arguments.
elemb : real(\*), in
Second argument in dot product.
**Details**
The interface of the accelerate package to the BLAS routine
zdotc appears for various packages, among them scipy as documented
in https://github.com/scipy/scipy/issues/1847. Using this workaround,
we loose possibly parallel implementations within the dot product.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function macdot_real(dim, elema, elemb) result(res)
integer, intent(in) :: dim
real(KIND=rKind), dimension(:), intent(in) :: elema
real(KIND=rKind), dimension(:), intent(in) :: elemb
real(KIND=rKind) :: res
res = sum((elema(:dim)) * elemb(:dim))
end function macdot_real
"""
return
[docs]def macdot_complex():
"""
fortran-function - July 2018 (dj)
Dot product implemented to avoid faulty ZDOTC implementation of Mac OS X
accelerate package.
**Arguments**
dim : INTEGER, in
Dimension of the dot product.
elema : complex(\*), in
First argument in dot product. Complex-conjugate taken
for complex-valued arguments.
elemb : complex(\*), in
Second argument in dot product.
**Details**
The interface of the accelerate package to the BLAS routine
zdotc appears for various packages, among them scipy as documented
in https://github.com/scipy/scipy/issues/1847. Using this workaround,
we loose possibly parallel implementations within the dot product.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function macdot_complex(dim, elema, elemb) result(res)
integer, intent(in) :: dim
complex(KIND=rKind), dimension(:), intent(in) :: elema
complex(KIND=rKind), dimension(:), intent(in) :: elemb
complex(KIND=rKind) :: res
res = sum(conjg(elema(:dim)) * elemb(:dim))
end function macdot_complex
"""
return
[docs]def equals_tensor_tensor():
"""
fortran-function - August 2017 (dj)
Compare if two tensors are equal.
**Arguments**
Tensa : TYPE(tensor), in
First tensor for comparison.
Tensb : TYPE(tensor), in
Second tensor for comparison.
tol : REAL, OPTIONAL, in
Tolerance for testing the actual entries of the tensor.
Default is 1e-14.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function equals_tensor_tensor(Tensa, Tensb, tol, &
errst) result(eq)
type(tensor), intent(in) :: Tensa
type(tensor), intent(in) :: Tensb
real(KIND=rKind), intent(in), optional :: tol
integer, intent(out), optional :: errst
logical :: eq
! Local variables
! ---------------
! error between the two arrays
real(KIND=rKind) :: err
!if(present(errst)) errst = 0
if(Tensa%rank /= Tensb%rank) then
eq = .false.
return
end if
if(any(Tensa%dl /= Tensb%dl)) then
eq = .false.
return
end if
err = maxval(abs(Tensa%elem(:Tensa%dim) - Tensb%elem(:Tensb%dim)))
if(present(tol)) then
eq = err < tol
else
eq = err < 1e-14_rKind
end if
end function equals_tensor_tensor
"""
return
[docs]def equals_tensorc_tensor():
"""
fortran-function - August 2017 (dj)
Compare if two tensors are equal.
**Arguments**
Tensa : TYPE(tensorc), in
First tensor for comparison.
Tensb : TYPE(tensor), in
Second tensor for comparison.
tol : REAL, OPTIONAL, in
Tolerance for testing the actual entries of the tensor.
Default is 1e-14.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function equals_tensorc_tensor(Tensa, Tensb, tol, &
errst) result(eq)
type(tensorc), intent(in) :: Tensa
type(tensor), intent(in) :: Tensb
real(KIND=rKind), intent(in), optional :: tol
integer, intent(out), optional :: errst
logical :: eq
! Local variables
! ---------------
! error between the two arrays
real(KIND=rKind) :: err
!if(present(errst)) errst = 0
if(Tensa%rank /= Tensb%rank) then
eq = .false.
return
end if
if(any(Tensa%dl /= Tensb%dl)) then
eq = .false.
return
end if
err = maxval(abs(Tensa%elem(:Tensa%dim) - Tensb%elem(:Tensb%dim)))
if(present(tol)) then
eq = err < tol
else
eq = err < 1e-14_rKind
end if
end function equals_tensorc_tensor
"""
return
[docs]def equals_tensor_tensorc():
"""
fortran-function - August 2017 (dj)
Compare if two tensors are equal.
**Arguments**
Tensa : TYPE(tensor), in
First tensor for comparison.
Tensb : TYPE(tensorc), in
Second tensor for comparison.
tol : REAL, OPTIONAL, in
Tolerance for testing the actual entries of the tensor.
Default is 1e-14.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function equals_tensor_tensorc(Tensa, Tensb, tol, &
errst) result(eq)
type(tensor), intent(in) :: Tensa
type(tensorc), intent(in) :: Tensb
real(KIND=rKind), intent(in), optional :: tol
integer, intent(out), optional :: errst
logical :: eq
! Local variables
! ---------------
! error between the two arrays
real(KIND=rKind) :: err
!if(present(errst)) errst = 0
if(Tensa%rank /= Tensb%rank) then
eq = .false.
return
end if
if(any(Tensa%dl /= Tensb%dl)) then
eq = .false.
return
end if
err = maxval(abs(Tensa%elem(:Tensa%dim) - Tensb%elem(:Tensb%dim)))
if(present(tol)) then
eq = err < tol
else
eq = err < 1e-14_rKind
end if
end function equals_tensor_tensorc
"""
return
[docs]def equals_tensorc_tensorc():
"""
fortran-function - August 2017 (dj)
Compare if two tensors are equal.
**Arguments**
Tensa : TYPE(tensorc), in
First tensor for comparison.
Tensb : TYPE(tensorc), in
Second tensor for comparison.
tol : REAL, OPTIONAL, in
Tolerance for testing the actual entries of the tensor.
Default is 1e-14.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function equals_tensorc_tensorc(Tensa, Tensb, tol, &
errst) result(eq)
type(tensorc), intent(in) :: Tensa
type(tensorc), intent(in) :: Tensb
real(KIND=rKind), intent(in), optional :: tol
integer, intent(out), optional :: errst
logical :: eq
! Local variables
! ---------------
! error between the two arrays
real(KIND=rKind) :: err
!if(present(errst)) errst = 0
if(Tensa%rank /= Tensb%rank) then
eq = .false.
return
end if
if(any(Tensa%dl /= Tensb%dl)) then
eq = .false.
return
end if
err = maxval(abs(Tensa%elem(:Tensa%dim) - Tensb%elem(:Tensb%dim)))
if(present(tol)) then
eq = err < tol
else
eq = err < 1e-14_rKind
end if
end function equals_tensorc_tensorc
"""
return
[docs]def expm_real_tensor():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : real, in
Additional scalar inside exp-function.
Tens : TYPE(tensor), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expm_real_tensor(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
real(KIND=rKind), intent(in) :: sc
type(tensor), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_real_tensor : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1088', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expm_real_real(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expm_real_tensor
"""
return
[docs]def expm_complex_tensor():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : complex, in
Additional scalar inside exp-function.
Tens : TYPE(tensor), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expm_complex_tensor(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
complex(KIND=rKind), intent(in) :: sc
type(tensor), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_complex_tensor : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1088', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expm_complex_real(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expm_complex_tensor
"""
return
[docs]def expm_real_tensorc():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : real, in
Additional scalar inside exp-function.
Tens : TYPE(tensorc), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expm_real_tensorc(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
real(KIND=rKind), intent(in) :: sc
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_real_tensorc : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1088', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expm_real_complex(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expm_real_tensorc
"""
return
[docs]def expm_complex_tensorc():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : complex, in
Additional scalar inside exp-function.
Tens : TYPE(tensorc), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expm_complex_tensorc(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
complex(KIND=rKind), intent(in) :: sc
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_complex_tensorc : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1088', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expm_complex_complex(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expm_complex_tensorc
"""
return
[docs]def expm_tensor_real_tensor():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying matrix.
This is a dummy.
**Arguments**
Texp : TYPE(tensor), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : real, in
Additional scalar inside exp-function.
Tens : TYPE(tensor), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Details**
The matrix exponential of a non-hermitian real matrix
can be complex and is therefore ill-defined. This methods just
raises and error.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expm_tensor_real_tensor(Texp, sc, Tens, &
last_row_idx, errst)
type(tensor), intent(out) :: Texp
real(KIND=rKind), intent(in) :: sc
type(tensor), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
errst = raise_error('exp_tensor_real_tensor'//&
': Non hermitian exp does not stay real.', &
99, 'Tensors_include.f90:1159', errst=errst)
end subroutine expm_tensor_real_tensor
"""
return
[docs]def expmh_real_tensor():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying hermitian
matrix.
**Arguments**
Texp : TYPE(tensor), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : real, in
Additional scalar inside exp-function.
Tens : TYPE(tensor), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Details**
Subroutines does not check if matrix is actually hermitian.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expmh_real_tensor(Texp, sc, Tens, last_row_idx, &
errst)
type(tensor), intent(out) :: Texp
real(KIND=rKind), intent(in) :: sc
type(tensor), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_real_tensor : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1232', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expmh_real_real(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expmh_real_tensor
"""
return
[docs]def expmh_complex_tensor():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying hermitian
matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : complex, in
Additional scalar inside exp-function.
Tens : TYPE(tensor), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Details**
Subroutines does not check if matrix is actually hermitian.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expmh_complex_tensor(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
complex(KIND=rKind), intent(in) :: sc
type(tensor), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_complex_tensor : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1232', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expmh_complex_real(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expmh_complex_tensor
"""
return
[docs]def expmh_real_tensorc():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying hermitian
matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : real, in
Additional scalar inside exp-function.
Tens : TYPE(tensorc), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Details**
Subroutines does not check if matrix is actually hermitian.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expmh_real_tensorc(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
real(KIND=rKind), intent(in) :: sc
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_real_tensorc : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1232', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expmh_real_complex(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expmh_real_tensorc
"""
return
[docs]def expmh_complex_tensorc():
"""
fortran-subroutine - August 2017 (dj)
Take the exponential of a tensor assuming an underlying hermitian
matrix.
**Arguments**
Texp : TYPE(tensorc), out
This is on exit the exponetial. Has rank and dimensions of
the input tensor.
sc : complex, in
Additional scalar inside exp-function.
Tens : TYPE(tensorc), inout
Take the exponential of this tensor.
last_row_idx : INTEGER, in
The dimension for the rows are calculated as
product(Tens%dl(:last_row_idx)). The remaining indices
build the column.
**Details**
Subroutines does not check if matrix is actually hermitian.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine expmh_complex_tensorc(Texp, sc, Tens, last_row_idx, &
errst)
type(tensorc), intent(out) :: Texp
complex(KIND=rKind), intent(in) :: sc
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: last_row_idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! matrix dimensions
integer :: d1, d2
d1 = product(Tens%dl(:last_row_idx))
d2 = product(Tens%dl(last_row_idx + 1:))
!if(d1 /= d2) then
! errst = raise_error('exp_complex_tensorc : '//&
! 'dimension mismatch.', 99, &
! 'Tensors_include.f90:1232', errst=errst)
! return
!end if
call create(Texp, Tens%dl)
call expmh_complex_complex(Texp%elem, sc, Tens%elem(:Tens%dim), &
d1, errst=errst)
end subroutine expmh_complex_tensorc
"""
return
[docs]def fuse_tensor():
"""
fortran-subroutine - October 2016 (dj)
Fuse indices of a tensor together to a single index.
**Arguments**
Tens : TYPE(tensor), inout
Merge certain indices of this tensor together to one index.
ind : INTEGER(\*), inout
Contains the positions of the indices to be merged together. Indices
must be ascending at this point.
**Details**
If indices are not neighbors, they are merged into the first
position.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine fuse_tensor(Tens, ind, errst)
type(tensor), intent(inout) :: Tens
integer, dimension(:), intent(in) :: ind
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping / indexing
integer :: ii, jj, kk
! old rank
integer :: orank
! copy of dimensions / permutations
integer, dimension(:), pointer :: darr
! Permutation to group merged indices together
integer, dimension(Tens%rank) :: perm
!if(present(errst)) errst = 0
!if(maxval(ind) > Tens%rank) then
! errst = raise_error('fuse_tensor: index exceeds '// &
! 'rank.', 2, 'Tensors_include.f90:1299', errst=errst)
! return
!end if
!if(any(ind(2:) - ind(:size(ind) - 1) < 0)) then
! errst = raise_error('fuse_tensor: indices not '// &
! 'ascending.', 5, 'Tensors_include.f90:1305', errst=errst)
! return
!end if
! Permute first to avoid unravelling this
call permute(Tens, errst=errst)
!if(prop_error('fuse_tensor: permute (1) failed,', &
! 'Tensors_include.f90:1312', errst=errst)) return
perm = [(ii, ii = 1, Tens%rank)]
do ii = 2, size(ind)
perm(ind(1) + ii - 1) = ind(ii)
end do
kk = ind(1) + size(ind)
do ii = 1, size(ind) - 1
do jj = ind(ii) + 1, ind(ii + 1) - 1
perm(kk) = jj
kk = kk + 1
end do
end do
do jj = ind(size(ind)) + 1, Tens%rank
perm(kk) = jj
kk = kk + 1
end do
Tens%idx = perm
call permute(Tens, errst=errst)
!if(prop_error('fuse_tensor: permute (2) failed,', &
! 'Tensors_include.f90:1335', errst=errst)) return
allocate(darr(Tens%rank))
darr = Tens%dl
orank = Tens%rank
Tens%rank = Tens%rank - size(ind) + 1
deallocate(Tens%dl, Tens%idx, Tens%mdl)
allocate(Tens%dl(Tens%rank), Tens%idx(Tens%rank), Tens%mdl(Tens%rank))
Tens%idx = [(ii, ii = 1, Tens%rank)]
! Reset dimensions
Tens%dl(1:ind(1) - 1) = darr(1:ind(1) - 1)
Tens%dl(ind(1)) = product(darr(ind(1):ind(1) + size(ind) - 1))
Tens%dl(ind(1) + 1:Tens%rank) = darr(ind(1) + size(ind):orank)
! Temporary solution, not really true
Tens%mdl = Tens%dl
deallocate(darr)
end subroutine fuse_tensor
"""
return
[docs]def fuse_tensorc():
"""
fortran-subroutine - October 2016 (dj)
Fuse indices of a tensor together to a single index.
**Arguments**
Tens : TYPE(tensorc), inout
Merge certain indices of this tensor together to one index.
ind : INTEGER(\*), inout
Contains the positions of the indices to be merged together. Indices
must be ascending at this point.
**Details**
If indices are not neighbors, they are merged into the first
position.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine fuse_tensorc(Tens, ind, errst)
type(tensorc), intent(inout) :: Tens
integer, dimension(:), intent(in) :: ind
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping / indexing
integer :: ii, jj, kk
! old rank
integer :: orank
! copy of dimensions / permutations
integer, dimension(:), pointer :: darr
! Permutation to group merged indices together
integer, dimension(Tens%rank) :: perm
!if(present(errst)) errst = 0
!if(maxval(ind) > Tens%rank) then
! errst = raise_error('fuse_tensorc: index exceeds '// &
! 'rank.', 2, 'Tensors_include.f90:1299', errst=errst)
! return
!end if
!if(any(ind(2:) - ind(:size(ind) - 1) < 0)) then
! errst = raise_error('fuse_tensorc: indices not '// &
! 'ascending.', 5, 'Tensors_include.f90:1305', errst=errst)
! return
!end if
! Permute first to avoid unravelling this
call permute(Tens, errst=errst)
!if(prop_error('fuse_tensorc: permute (1) failed,', &
! 'Tensors_include.f90:1312', errst=errst)) return
perm = [(ii, ii = 1, Tens%rank)]
do ii = 2, size(ind)
perm(ind(1) + ii - 1) = ind(ii)
end do
kk = ind(1) + size(ind)
do ii = 1, size(ind) - 1
do jj = ind(ii) + 1, ind(ii + 1) - 1
perm(kk) = jj
kk = kk + 1
end do
end do
do jj = ind(size(ind)) + 1, Tens%rank
perm(kk) = jj
kk = kk + 1
end do
Tens%idx = perm
call permute(Tens, errst=errst)
!if(prop_error('fuse_tensorc: permute (2) failed,', &
! 'Tensors_include.f90:1335', errst=errst)) return
allocate(darr(Tens%rank))
darr = Tens%dl
orank = Tens%rank
Tens%rank = Tens%rank - size(ind) + 1
deallocate(Tens%dl, Tens%idx, Tens%mdl)
allocate(Tens%dl(Tens%rank), Tens%idx(Tens%rank), Tens%mdl(Tens%rank))
Tens%idx = [(ii, ii = 1, Tens%rank)]
! Reset dimensions
Tens%dl(1:ind(1) - 1) = darr(1:ind(1) - 1)
Tens%dl(ind(1)) = product(darr(ind(1):ind(1) + size(ind) - 1))
Tens%dl(ind(1) + 1:Tens%rank) = darr(ind(1) + size(ind):orank)
! Temporary solution, not really true
Tens%mdl = Tens%dl
deallocate(darr)
end subroutine fuse_tensorc
"""
return
[docs]def fuse_all_tensor():
"""
fortran-subroutine - October 2017 (dj)
Fuse equally sized subsets of links containing in total all links.
**Arguments**
Tens : TYPE(tensor), inout
Fuse links in this tensor.
idx : INTEGER(\*, \*), in
Each column contain the links to be fused. The number of rows times
the number of columns must correspond to the rank of the tensor.
ii-th column is ii-th new index.
method : CHARACTER, in
Dummy flag to comply with the interface for the symmetric tensors.
ordered : LOGICAL, OPTIONAL, in
If flag is true, the permutation of the indices does not have to
be executed. Default to false (permute links).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine fuse_all_tensor(Tens, idx, method, ordered, errst)
type(tensor), intent(inout) :: Tens
integer, dimension(:, :), intent(in) :: idx
character, intent(in) :: method
logical, intent(in), optional :: ordered
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! for indexing while looping
integer :: i1, i2
! duplette for optional argument
logical :: ordered_
! permutation to put fused indices next to each other
integer, dimension(:), allocatable :: perm
! new dimension of the tensor
integer, dimension(:), allocatable :: newdl
! Avoid unused variable message for intentionally unused variables
!if(.false.) print *, method
!if(present(errst)) errst = 0
!if(size(idx) /= Tens%rank) then
! errst = raise_error('fuse_all_tensor : rank mismatch.', &
! 99, 'Tensors_include.f90:1430', errst=errst)
! return
!end if
ordered_ = .false.
if(present(ordered)) ordered_ = ordered
if(.not. ordered_) then
allocate(perm(size(idx)))
perm = reshape(idx, [size(idx)])
call transposed(Tens, perm, doperm=.true.)
deallocate(perm)
end if
Tens%rank = size(idx, 2)
allocate(newdl(Tens%rank))
i1 = 1
i2 = size(idx, 1)
do ii = 1, Tens%rank
newdl(ii) = product(Tens%dl(i1:i2))
i1 = i1 + size(idx, 1)
i2 = i2 + size(idx, 1)
end do
deallocate(Tens%dl, Tens%idx)
allocate(Tens%dl(Tens%rank), Tens%idx(Tens%rank))
Tens%dl = newdl
Tens%idx = [(ii, ii = 1, Tens%rank)]
! Figure out maximal dimension
i1 = 1
i2 = size(idx, 1)
do ii = 1, Tens%rank
newdl(ii) = product(Tens%mdl(i1:i2))
i1 = i1 + size(idx, 1)
i2 = i2 + size(idx, 1)
end do
deallocate(Tens%mdl)
allocate(Tens%mdl(Tens%rank))
Tens%mdl = newdl
deallocate(newdl)
end subroutine fuse_all_tensor
"""
return
[docs]def fuse_all_tensorc():
"""
fortran-subroutine - October 2017 (dj)
Fuse equally sized subsets of links containing in total all links.
**Arguments**
Tens : TYPE(tensorc), inout
Fuse links in this tensor.
idx : INTEGER(\*, \*), in
Each column contain the links to be fused. The number of rows times
the number of columns must correspond to the rank of the tensor.
ii-th column is ii-th new index.
method : CHARACTER, in
Dummy flag to comply with the interface for the symmetric tensors.
ordered : LOGICAL, OPTIONAL, in
If flag is true, the permutation of the indices does not have to
be executed. Default to false (permute links).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine fuse_all_tensorc(Tens, idx, method, ordered, errst)
type(tensorc), intent(inout) :: Tens
integer, dimension(:, :), intent(in) :: idx
character, intent(in) :: method
logical, intent(in), optional :: ordered
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! for indexing while looping
integer :: i1, i2
! duplette for optional argument
logical :: ordered_
! permutation to put fused indices next to each other
integer, dimension(:), allocatable :: perm
! new dimension of the tensor
integer, dimension(:), allocatable :: newdl
! Avoid unused variable message for intentionally unused variables
!if(.false.) print *, method
!if(present(errst)) errst = 0
!if(size(idx) /= Tens%rank) then
! errst = raise_error('fuse_all_tensorc : rank mismatch.', &
! 99, 'Tensors_include.f90:1430', errst=errst)
! return
!end if
ordered_ = .false.
if(present(ordered)) ordered_ = ordered
if(.not. ordered_) then
allocate(perm(size(idx)))
perm = reshape(idx, [size(idx)])
call transposed(Tens, perm, doperm=.true.)
deallocate(perm)
end if
Tens%rank = size(idx, 2)
allocate(newdl(Tens%rank))
i1 = 1
i2 = size(idx, 1)
do ii = 1, Tens%rank
newdl(ii) = product(Tens%dl(i1:i2))
i1 = i1 + size(idx, 1)
i2 = i2 + size(idx, 1)
end do
deallocate(Tens%dl, Tens%idx)
allocate(Tens%dl(Tens%rank), Tens%idx(Tens%rank))
Tens%dl = newdl
Tens%idx = [(ii, ii = 1, Tens%rank)]
! Figure out maximal dimension
i1 = 1
i2 = size(idx, 1)
do ii = 1, Tens%rank
newdl(ii) = product(Tens%mdl(i1:i2))
i1 = i1 + size(idx, 1)
i2 = i2 + size(idx, 1)
end do
deallocate(Tens%mdl)
allocate(Tens%mdl(Tens%rank))
Tens%mdl = newdl
deallocate(newdl)
end subroutine fuse_all_tensorc
"""
return
[docs]def gaxpy_tensor_real_tensor():
"""
fortran-subroutine - March 2016 (updated dj)
Y = Y + a * X for real a
**Arguments**
Yy : TYPE(tensor), inout
On entry, tensor. On exit tensor times scalar added to this tensor,
sc : real, in
Scalar multiplied with tensor X.
Xx : TYPE(tensor), inout
Tensor added to tensor Y. Is scaled with constant `sc`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine gaxpy_tensor_real_tensor(Yy, sc, Xx, errst)
type(tensor), intent(inout) :: Yy
real(KIND=rKind), intent(in) :: sc
type(tensor), intent(in) :: Xx
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(Yy%rank /= Xx%rank) then
! errst = raise_error('gaxpy_tensor_real_'// &
! 'tensor: rank mismatch.', 2, &
! 'Tensors_include.f90:1531', errst=errst)
! return
!end if
!if((Yy%rank /= 0) .and. any(Yy%dl /= Xx%dl)) then
! errst = raise_error('gaxpy_tensor_real_'// &
! 'tensor: dimension mismatch.', 3, &
! 'Tensors_include.f90:1538', errst=errst)
! return
!end if
call DAXPY(Yy%dim, sc, Xx%elem, 1, Yy%elem, 1)
!Yy%elem(:Yy%dim) = Yy%elem(:Yy%dim) + sc * Xx%elem(:Yy%dim)
end subroutine gaxpy_tensor_real_tensor
"""
return
[docs]def gaxpy_tensorc_complex_tensorc():
"""
fortran-subroutine - March 2016 (updated dj)
Y = Y + a * X for complex a
**Arguments**
Yy : TYPE(tensorc), inout
On entry, tensor. On exit tensor times scalar added to this tensor,
sc : complex, in
Scalar multiplied with tensor X.
Xx : TYPE(tensorc), inout
Tensor added to tensor Y. Is scaled with constant `sc`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine gaxpy_tensorc_complex_tensorc(Yy, sc, Xx, errst)
type(tensorc), intent(inout) :: Yy
complex(KIND=rKind), intent(in) :: sc
type(tensorc), intent(in) :: Xx
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(Yy%rank /= Xx%rank) then
! errst = raise_error('gaxpy_tensorc_complex_'// &
! 'tensorc: rank mismatch.', 2, &
! 'Tensors_include.f90:1531', errst=errst)
! return
!end if
!if((Yy%rank /= 0) .and. any(Yy%dl /= Xx%dl)) then
! errst = raise_error('gaxpy_tensorc_complex_'// &
! 'tensorc: dimension mismatch.', 3, &
! 'Tensors_include.f90:1538', errst=errst)
! return
!end if
call ZAXPY(Yy%dim, sc, Xx%elem, 1, Yy%elem, 1)
!Yy%elem(:Yy%dim) = Yy%elem(:Yy%dim) + sc * Xx%elem(:Yy%dim)
end subroutine gaxpy_tensorc_complex_tensorc
"""
return
[docs]def gaxpy_tensorc_real_tensor():
"""
fortran-subroutine - March 2016 (updated dj)
Y = Y + a * X for real a
**Arguments**
Yy : TYPE(tensorc), inout
On entry, tensor. On exit tensor times scalar added to this tensor,
sc : real, in
Scalar multiplied with tensor X.
Xx : TYPE(tensor), inout
Tensor added to tensor Y. Is scaled with constant `sc`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine gaxpy_tensorc_real_tensor(Yy, sc, Xx, errst)
type(tensorc), intent(inout) :: Yy
real(KIND=rKind), intent(in) :: sc
type(tensor), intent(in) :: Xx
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(Yy%rank /= Xx%rank) then
! errst = raise_error('gaxpy_tensorc_real_'// &
! 'tensor: rank mismatch.', 2, &
! 'Tensors_include.f90:1531', errst=errst)
! return
!end if
!if((Yy%rank /= 0) .and. any(Yy%dl /= Xx%dl)) then
! errst = raise_error('gaxpy_tensorc_real_'// &
! 'tensor: dimension mismatch.', 3, &
! 'Tensors_include.f90:1538', errst=errst)
! return
!end if
! (Yy%dim, sc, Xx%elem, 1, Yy%elem, 1)
Yy%elem(:Yy%dim) = Yy%elem(:Yy%dim) + sc * Xx%elem(:Yy%dim)
end subroutine gaxpy_tensorc_real_tensor
"""
return
[docs]def gaxpy_tensorc_complex_tensor():
"""
fortran-subroutine - March 2016 (updated dj)
Y = Y + a * X for complex a
**Arguments**
Yy : TYPE(tensorc), inout
On entry, tensor. On exit tensor times scalar added to this tensor,
sc : complex, in
Scalar multiplied with tensor X.
Xx : TYPE(tensor), inout
Tensor added to tensor Y. Is scaled with constant `sc`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine gaxpy_tensorc_complex_tensor(Yy, sc, Xx, errst)
type(tensorc), intent(inout) :: Yy
complex(KIND=rKind), intent(in) :: sc
type(tensor), intent(in) :: Xx
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(Yy%rank /= Xx%rank) then
! errst = raise_error('gaxpy_tensorc_complex_'// &
! 'tensor: rank mismatch.', 2, &
! 'Tensors_include.f90:1531', errst=errst)
! return
!end if
!if((Yy%rank /= 0) .and. any(Yy%dl /= Xx%dl)) then
! errst = raise_error('gaxpy_tensorc_complex_'// &
! 'tensor: dimension mismatch.', 3, &
! 'Tensors_include.f90:1538', errst=errst)
! return
!end if
! (Yy%dim, sc, Xx%elem, 1, Yy%elem, 1)
Yy%elem(:Yy%dim) = Yy%elem(:Yy%dim) + sc * Xx%elem(:Yy%dim)
end subroutine gaxpy_tensorc_complex_tensor
"""
return
[docs]def gaxpy_tensorc_real_tensorc():
"""
fortran-subroutine - March 2016 (updated dj)
Y = Y + a * X for real a
**Arguments**
Yy : TYPE(tensorc), inout
On entry, tensor. On exit tensor times scalar added to this tensor,
sc : real, in
Scalar multiplied with tensor X.
Xx : TYPE(tensorc), inout
Tensor added to tensor Y. Is scaled with constant `sc`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine gaxpy_tensorc_real_tensorc(Yy, sc, Xx, errst)
type(tensorc), intent(inout) :: Yy
real(KIND=rKind), intent(in) :: sc
type(tensorc), intent(in) :: Xx
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
!if(Yy%rank /= Xx%rank) then
! errst = raise_error('gaxpy_tensorc_real_'// &
! 'tensorc: rank mismatch.', 2, &
! 'Tensors_include.f90:1531', errst=errst)
! return
!end if
!if((Yy%rank /= 0) .and. any(Yy%dl /= Xx%dl)) then
! errst = raise_error('gaxpy_tensorc_real_'// &
! 'tensorc: dimension mismatch.', 3, &
! 'Tensors_include.f90:1538', errst=errst)
! return
!end if
! (Yy%dim, sc, Xx%elem, 1, Yy%elem, 1)
Yy%elem(:Yy%dim) = Yy%elem(:Yy%dim) + sc * Xx%elem(:Yy%dim)
end subroutine gaxpy_tensorc_real_tensorc
"""
return
[docs]def get_scalar_tensor():
"""
fortran-function - July 2017 (dj)
Get the scalar if the tensor represents a rank-0 tensor.
**Arguments**
Qtens : TYPE(Qtensor), in
Return the scalar entry of a rank-0 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function get_scalar_tensor(Tens, errst) result (sc)
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
real(KIND=rKind) :: sc
! No local variables
! ------------------
!if(present(errst)) errst = 0
sc = Tens%elem(1)
end function get_scalar_tensor
"""
return
[docs]def get_scalar_tensorc():
"""
fortran-function - July 2017 (dj)
Get the scalar if the tensor represents a rank-0 tensor.
**Arguments**
Qtens : TYPE(Qtensorc), in
Return the scalar entry of a rank-0 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function get_scalar_tensorc(Tens, errst) result (sc)
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
complex(KIND=rKind) :: sc
! No local variables
! ------------------
!if(present(errst)) errst = 0
sc = Tens%elem(1)
end function get_scalar_tensorc
"""
return
[docs]def has_nan_tensor():
"""
fortran-function - October 2017 (dj)
Check if the tensor has NAN in the actual array.
**Arguments**
Tens : TYPE(tensor), in
Check NAN on the entries of this tensor. Check does not
apply to rank, dimensions, etc.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function has_nan_tensor(Tens, errst) result(hasnan)
use ieee_arithmetic, only : ieee_is_nan
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
logical :: hasnan
! No local variables
! ------------------
!if(present(errst)) errst = 0
hasnan = ieee_is_nan(sum(abs(Tens%elem(:Tens%dim))))
end function has_nan_tensor
"""
return
[docs]def has_nan_tensorc():
"""
fortran-function - October 2017 (dj)
Check if the tensor has NAN in the actual array.
**Arguments**
Tens : TYPE(tensorc), in
Check NAN on the entries of this tensor. Check does not
apply to rank, dimensions, etc.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function has_nan_tensorc(Tens, errst) result(hasnan)
use ieee_arithmetic, only : ieee_is_nan
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
logical :: hasnan
! No local variables
! ------------------
!if(present(errst)) errst = 0
hasnan = ieee_is_nan(sum(abs(Tens%elem(:Tens%dim))))
end function has_nan_tensorc
"""
return
[docs]def is_eye_tensor():
"""
fortran-function - August 2017 (dj)
Check if the tensor is an identity.
**Arguments**
Tens : TYPE(tensor), in
This tensor is checked if it is an identity.
**Details**
We define the identity for any rank as a tensor
with ones for i_1 = i_2 = ... = i_n, and otherwise zeros.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function is_eye_tensor(Tens, errst) result(res)
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
logical :: res
! Local variables
! ---------------
! for looping
integer :: ii
! stride for to get next element
integer :: stride
! index running over diagonal elements
integer :: idx
! Copy for convenience
type(tensor) :: Tmp
call copy(Tmp, Tens)
res = .true.
stride = 1
do ii = 1, (Tens%rank - 1)
stride = stride + product(Tens%dl(:ii))
end do
idx = 1
do ii = 1, minval(Tens%dl)
if(abs(Tmp%elem(idx) - 1.0_rKind) > 1e-14_rKind) then
res = .false.
call destroy(Tmp)
return
else
Tmp%elem(idx) = 0.0_rKind
end if
idx = idx + stride
end do
if(maxval(abs(Tmp%elem(:Tmp%dim))) > 1e-14) then
res = .false.
end if
call destroy(Tmp)
end function is_eye_tensor
"""
return
[docs]def is_eye_tensorc():
"""
fortran-function - August 2017 (dj)
Check if the tensor is an identity.
**Arguments**
Tens : TYPE(tensorc), in
This tensor is checked if it is an identity.
**Details**
We define the identity for any rank as a tensor
with ones for i_1 = i_2 = ... = i_n, and otherwise zeros.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function is_eye_tensorc(Tens, errst) result(res)
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
logical :: res
! Local variables
! ---------------
! for looping
integer :: ii
! stride for to get next element
integer :: stride
! index running over diagonal elements
integer :: idx
! Copy for convenience
type(tensorc) :: Tmp
call copy(Tmp, Tens)
res = .true.
stride = 1
do ii = 1, (Tens%rank - 1)
stride = stride + product(Tens%dl(:ii))
end do
idx = 1
do ii = 1, minval(Tens%dl)
if(abs(Tmp%elem(idx) - 1.0_rKind) > 1e-14_rKind) then
res = .false.
call destroy(Tmp)
return
else
Tmp%elem(idx) = 0.0_rKind
end if
idx = idx + stride
end do
if(maxval(abs(Tmp%elem(:Tmp%dim))) > 1e-14) then
res = .false.
end if
call destroy(Tmp)
end function is_eye_tensorc
"""
return
[docs]def is_set_tensor():
"""
fortran-function - August 2017 (dj)
Check if a tensor is in use.
**Arguments**
Tens : TYPE(tensor), in
Check if the tensor was used, i.e. a call to create. The check
is on the rank of the tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function is_set_tensor(Tens, errst) result(res)
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
logical :: res
! No local variables
! ------------------
!if(present(errst)) errst = 0
res = .not. (Tens%rank == -1)
end function is_set_tensor
"""
return
[docs]def is_set_tensorc():
"""
fortran-function - August 2017 (dj)
Check if a tensor is in use.
**Arguments**
Tens : TYPE(tensorc), in
Check if the tensor was used, i.e. a call to create. The check
is on the rank of the tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function is_set_tensorc(Tens, errst) result(res)
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
logical :: res
! No local variables
! ------------------
!if(present(errst)) errst = 0
res = .not. (Tens%rank == -1)
end function is_set_tensorc
"""
return
[docs]def kron_tensor_tensor():
"""
fortran-subroutine - May 2017 (dj)
Kronecker product bipartions on two tensors.
**Arguments**
Tout : TYPE(tensor), inout
Storing the result of the Kronecker product.
Tl : TYPE(tensor), inout
Left tensor in the Kronecker product.
Tr : TYPE(tensor), inout
Right tensor in the Kronecker product.
llr : INTEGER, in
Last index contributing to the rows, all over 1:llr. This is
for the left tensor Tl.
lrr : INTEGER, in
Last index contributing to the rows, all over 1:lrr. This is
for the left tensor Tr.
transl : CHARACTER, in
Defines the transformation executed on the left matrix :math:`Tl`
when taking the tensor product. Either 'N' (none), 'T' (transposed),
or 'D' (dagger / conjugate transposed).
transr : CHARACTER, in
Defines the transformation executed on the left matrix :math:`Tr`
when taking the tensor product. Either 'N' (none), 'T' (transposed),
or 'D' (dagger / conjugate transposed).
op : CHARACTER, in
Either 'N' for usual Kronecker product, or '+' for adding
to the result.
**Details**
The output is always a rank-4 tensor. The intended use is
for rank-2 tensors.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine kron_tensor_tensor(Tout, Tl, Tr, llr, lrr, &
transl, transr, op, errst)
type(tensor), intent(inout) :: Tout
type(tensor), intent(in) :: Tl, Tr
integer, intent(in) :: llr, lrr
character, intent(in) :: transl, transr, op
integer, intent(out), optional :: errst
! Local variables
! ---------------
! dimension for rows and columns of left and right matrix
integer :: rl, cl, rr, cr
!if(present(errst)) errst = 0
rl = product(Tl%dl(:llr))
cl = product(Tl%dl(llr+1:))
rr = product(Tr%dl(:lrr))
cr = product(Tr%dl(lrr+1:))
if(op == 'N') then
call create(Tout, [rl, rr, cl, cr])
end if
call kron_real_real(Tout%elem, Tl%elem, Tr%elem, &
rl, cl, rr, cr, transl, transr, op=op)
end subroutine kron_tensor_tensor
"""
return
[docs]def kron_tensorc_tensorc():
"""
fortran-subroutine - May 2017 (dj)
Kronecker product bipartions on two tensors.
**Arguments**
Tout : TYPE(tensorc), inout
Storing the result of the Kronecker product.
Tl : TYPE(tensorc), inout
Left tensor in the Kronecker product.
Tr : TYPE(tensorc), inout
Right tensor in the Kronecker product.
llr : INTEGER, in
Last index contributing to the rows, all over 1:llr. This is
for the left tensor Tl.
lrr : INTEGER, in
Last index contributing to the rows, all over 1:lrr. This is
for the left tensor Tr.
transl : CHARACTER, in
Defines the transformation executed on the left matrix :math:`Tl`
when taking the tensor product. Either 'N' (none), 'T' (transposed),
or 'D' (dagger / conjugate transposed).
transr : CHARACTER, in
Defines the transformation executed on the left matrix :math:`Tr`
when taking the tensor product. Either 'N' (none), 'T' (transposed),
or 'D' (dagger / conjugate transposed).
op : CHARACTER, in
Either 'N' for usual Kronecker product, or '+' for adding
to the result.
**Details**
The output is always a rank-4 tensor. The intended use is
for rank-2 tensors.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine kron_tensorc_tensorc(Tout, Tl, Tr, llr, lrr, &
transl, transr, op, errst)
type(tensorc), intent(inout) :: Tout
type(tensorc), intent(in) :: Tl, Tr
integer, intent(in) :: llr, lrr
character, intent(in) :: transl, transr, op
integer, intent(out), optional :: errst
! Local variables
! ---------------
! dimension for rows and columns of left and right matrix
integer :: rl, cl, rr, cr
!if(present(errst)) errst = 0
rl = product(Tl%dl(:llr))
cl = product(Tl%dl(llr+1:))
rr = product(Tr%dl(:lrr))
cr = product(Tr%dl(lrr+1:))
if(op == 'N') then
call create(Tout, [rl, rr, cl, cr])
end if
call kron_complex_complex(Tout%elem, Tl%elem, Tr%elem, &
rl, cl, rr, cr, transl, transr, op=op)
end subroutine kron_tensorc_tensorc
"""
return
[docs]def kron_tensorc_tensor():
"""
fortran-subroutine - May 2017 (dj)
Kronecker product bipartions on two tensors.
**Arguments**
Tout : TYPE(tensorc), inout
Storing the result of the Kronecker product.
Tl : TYPE(tensor), inout
Left tensor in the Kronecker product.
Tr : TYPE(tensor), inout
Right tensor in the Kronecker product.
llr : INTEGER, in
Last index contributing to the rows, all over 1:llr. This is
for the left tensor Tl.
lrr : INTEGER, in
Last index contributing to the rows, all over 1:lrr. This is
for the left tensor Tr.
transl : CHARACTER, in
Defines the transformation executed on the left matrix :math:`Tl`
when taking the tensor product. Either 'N' (none), 'T' (transposed),
or 'D' (dagger / conjugate transposed).
transr : CHARACTER, in
Defines the transformation executed on the left matrix :math:`Tr`
when taking the tensor product. Either 'N' (none), 'T' (transposed),
or 'D' (dagger / conjugate transposed).
op : CHARACTER, in
Either 'N' for usual Kronecker product, or '+' for adding
to the result.
**Details**
The output is always a rank-4 tensor. The intended use is
for rank-2 tensors.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine kron_tensorc_tensor(Tout, Tl, Tr, llr, lrr, &
transl, transr, op, errst)
type(tensorc), intent(inout) :: Tout
type(tensor), intent(in) :: Tl, Tr
integer, intent(in) :: llr, lrr
character, intent(in) :: transl, transr, op
integer, intent(out), optional :: errst
! Local variables
! ---------------
! dimension for rows and columns of left and right matrix
integer :: rl, cl, rr, cr
!if(present(errst)) errst = 0
rl = product(Tl%dl(:llr))
cl = product(Tl%dl(llr+1:))
rr = product(Tr%dl(:lrr))
cr = product(Tr%dl(lrr+1:))
if(op == 'N') then
call create(Tout, [rl, rr, cl, cr])
end if
call kron_complex_real(Tout%elem, Tl%elem, Tr%elem, &
rl, cl, rr, cr, transl, transr, op=op)
end subroutine kron_tensorc_tensor
"""
return
[docs]def maxlineardim_tensor():
"""
fortran-function - June 2017 (updated dj)
Return the maximal dimension of all n indices in the
rank-n tensor, or selected indices.
**Arguments**
Tens : TYPE(tensor), in
Get the maximal dimension of this tensor.
idx : INTEGER(\*), in
Search only over these indices of the rank-n tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function maxlineardim_tensor(Tens, idx) result (mld)
type(tensor), intent(in) :: Tens
integer, dimension(:), optional :: idx
integer :: mld
if(present(idx)) then
mld = maxval(Tens%dl(idx))
else
mld = maxval(Tens%dl)
end if
end function maxlineardim_tensor
"""
return
[docs]def maxlineardim_tensorc():
"""
fortran-function - June 2017 (updated dj)
Return the maximal dimension of all n indices in the
rank-n tensor, or selected indices.
**Arguments**
Tens : TYPE(tensorc), in
Get the maximal dimension of this tensor.
idx : INTEGER(\*), in
Search only over these indices of the rank-n tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function maxlineardim_tensorc(Tens, idx) result (mld)
type(tensorc), intent(in) :: Tens
integer, dimension(:), optional :: idx
integer :: mld
if(present(idx)) then
mld = maxval(Tens%dl(idx))
else
mld = maxval(Tens%dl)
end if
end function maxlineardim_tensorc
"""
return
[docs]def maxvalue_tensor():
"""
fortran-function - June 2017 (dj)
Return the maximal value of all entries in the tensor.
**Arguments**
Tens : TYPE(tensor), in
Get the maximal entry of this tensor. For complex tensors,
the absolute value is considered.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function maxvalue_tensor(Tens) result(val)
type(tensor), intent(in) :: Tens
real(KIND=rKind) :: val
! Local variables
! ---------------
! index of maximum
integer :: idx
idx = maxloc((Tens%elem(:Tens%dim)), 1)
val = Tens%elem(idx)
end function maxvalue_tensor
"""
return
[docs]def maxvalue_tensorc():
"""
fortran-function - June 2017 (dj)
Return the maximal value of all entries in the tensor.
**Arguments**
Tens : TYPE(tensorc), in
Get the maximal entry of this tensor. For complex tensors,
the absolute value is considered.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function maxvalue_tensorc(Tens) result(val)
type(tensorc), intent(in) :: Tens
complex(KIND=rKind) :: val
! Local variables
! ---------------
! index of maximum
integer :: idx
idx = maxloc(abs(Tens%elem(:Tens%dim)), 1)
val = Tens%elem(idx)
end function maxvalue_tensorc
"""
return
[docs]def norm_tensor():
"""
fortran-function - October 2016 (update dj)
Compute <A, A> for a tensor.
**Arguments**
Tens : TYPE(tensor), in
Get norm of this tensor <A, A>.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function norm_tensor(Tens, errst) result(norm)
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
real(KIND=rKind) :: norm
! External functions
! ------------------
! goes into DDOT,ZDOTC
real(KIND=rKind), external :: DDOT
!if(present(errst)) errst = 0
norm = real(DDOT(Tens%dim, Tens%elem, 1, &
Tens%elem, 1), KIND=rKind)
!norm = real(macdot(Tens%dim, Tens%elem, Tens%elem), &
! KIND=rKind)
end function norm_tensor
"""
return
[docs]def norm_tensorc():
"""
fortran-function - October 2016 (update dj)
Compute <A, A> for a tensor.
**Arguments**
Tens : TYPE(tensorc), in
Get norm of this tensor <A, A>.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function norm_tensorc(Tens, errst) result(norm)
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
real(KIND=rKind) :: norm
! External functions
! ------------------
! goes into DDOT,ZDOTC
complex(KIND=rKind), external :: ZDOTC
!if(present(errst)) errst = 0
norm = real(ZDOTC(Tens%dim, Tens%elem, 1, &
Tens%elem, 1), KIND=rKind)
!norm = real(macdot(Tens%dim, Tens%elem, Tens%elem), &
! KIND=rKind)
end function norm_tensorc
"""
return
[docs]def permute_tensor():
"""
fortran-subroutine - July 2016 (dj)
Bring the tensor in the actual order defined in idx in the memory.
It is "COME FROM" indexing.
**Arguments**
Tens : TYPE(tensor), inout
Permute the indices of this tensor into the order stored in idx.
**Details**
We divide the indices into three blocks (A, B, C), where the
first and the last block contains indices not changing their position.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine permute_tensor(Tens, errst)
type(tensor), intent(inout) :: Tens
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii, jj, kk
! index in new array
integer :: idx
! number of indices in each block (A, B, C)
integer :: na, nb, nc
! total dimension of indices in each block (A, B, C)
integer :: da, db, dc
! product of da * db
integer :: dab
! pos in old and new array
integer, dimension(:), pointer :: posold, newdim, cumdim
! temporary memory
real(KIND=rKIND), dimension(:), pointer :: pelem
!if(present(errst)) errst = 0
!if(size(Tens%idx, 1) > Tens%rank) then
! errst = raise_error('permute_tensor: dimension '// &
! 'mismatch.', 99, 'Tensors_include.f90:2005', errst=errst)
! return
!end if
! Find indices in block A
na = 0
da = 1
do ii = 1, Tens%rank
if(ii == Tens%idx(ii)) then
na = na + 1
da = da * Tens%dl(ii)
else
exit
end if
end do
! If block size is dim - no permutation
if(na == Tens%rank) return
allocate(pelem(size(Tens%elem, 1)))
! Find indices in block C
nc = 0
dc = 1
do ii = Tens%rank, 1, -1
if(ii == Tens%idx(ii)) then
nc = nc + 1
dc = dc * Tens%dl(ii)
else
exit
end if
end do
! Set block B
nb = Tens%rank - na - nc
db = 1
do ii = 1, nb
db = db * Tens%dl(na + ii)
end do
! Set dimension in block B after permutation
allocate(cumdim(nb), newdim(nb))
newdim(1) = da
do ii = na + 1, na + nb - 1
newdim(ii - na + 1) = newdim(ii - na) * Tens%dl(Tens%idx(ii))
end do
do ii = na + 1, na + nb
cumdim(Tens%idx(ii) - na) = newdim(ii - na)
end do
! Arrays for old position
allocate(posold(nb))
posold = 0
posold(1) = -1
! Loop over indices to be permuted
! --------------------------------
if((na /= 0) .and. (nc /= 0)) then
! Unchagend indices in front and back
dab = da * db
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
kk = 0
do jj = 1, dc
pelem(idx:idx + da - 1) = &
Tens%elem(kk + da * (ii - 1) + 1:kk + da * ii)
idx = idx + dab
kk = kk + dab
end do
end do
elseif(na /= 0) then
! Unchanged indices in front
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
pelem(idx:idx + da - 1) = Tens%elem(da * (ii - 1) + 1:da * ii)
end do
elseif(nc /= 0) then
! Unchanged indices in back
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
kk = ii
do jj = 1, dc
pelem(idx) = Tens%elem(kk)
idx = idx + db
kk = kk + db
end do
end do
else
! No unchanged indices in front or back
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
pelem(idx) = Tens%elem(ii)
end do
end if
deallocate(Tens%elem)
Tens%elem => pelem
pelem => null()
!Tens%elem(:Tens%dim) = pelem
Tens%dl = Tens%dl(Tens%idx)
Tens%idx = [(ii, ii = 1, Tens%rank)]
deallocate(newdim, posold, cumdim)
!deallocate(pelem, newdim, posold, cumdim)
end subroutine permute_tensor
"""
return
[docs]def permute_tensorc():
"""
fortran-subroutine - July 2016 (dj)
Bring the tensor in the actual order defined in idx in the memory.
It is "COME FROM" indexing.
**Arguments**
Tens : TYPE(tensorc), inout
Permute the indices of this tensor into the order stored in idx.
**Details**
We divide the indices into three blocks (A, B, C), where the
first and the last block contains indices not changing their position.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine permute_tensorc(Tens, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii, jj, kk
! index in new array
integer :: idx
! number of indices in each block (A, B, C)
integer :: na, nb, nc
! total dimension of indices in each block (A, B, C)
integer :: da, db, dc
! product of da * db
integer :: dab
! pos in old and new array
integer, dimension(:), pointer :: posold, newdim, cumdim
! temporary memory
complex(KIND=rKIND), dimension(:), pointer :: pelem
!if(present(errst)) errst = 0
!if(size(Tens%idx, 1) > Tens%rank) then
! errst = raise_error('permute_tensorc: dimension '// &
! 'mismatch.', 99, 'Tensors_include.f90:2005', errst=errst)
! return
!end if
! Find indices in block A
na = 0
da = 1
do ii = 1, Tens%rank
if(ii == Tens%idx(ii)) then
na = na + 1
da = da * Tens%dl(ii)
else
exit
end if
end do
! If block size is dim - no permutation
if(na == Tens%rank) return
allocate(pelem(size(Tens%elem, 1)))
! Find indices in block C
nc = 0
dc = 1
do ii = Tens%rank, 1, -1
if(ii == Tens%idx(ii)) then
nc = nc + 1
dc = dc * Tens%dl(ii)
else
exit
end if
end do
! Set block B
nb = Tens%rank - na - nc
db = 1
do ii = 1, nb
db = db * Tens%dl(na + ii)
end do
! Set dimension in block B after permutation
allocate(cumdim(nb), newdim(nb))
newdim(1) = da
do ii = na + 1, na + nb - 1
newdim(ii - na + 1) = newdim(ii - na) * Tens%dl(Tens%idx(ii))
end do
do ii = na + 1, na + nb
cumdim(Tens%idx(ii) - na) = newdim(ii - na)
end do
! Arrays for old position
allocate(posold(nb))
posold = 0
posold(1) = -1
! Loop over indices to be permuted
! --------------------------------
if((na /= 0) .and. (nc /= 0)) then
! Unchagend indices in front and back
dab = da * db
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
kk = 0
do jj = 1, dc
pelem(idx:idx + da - 1) = &
Tens%elem(kk + da * (ii - 1) + 1:kk + da * ii)
idx = idx + dab
kk = kk + dab
end do
end do
elseif(na /= 0) then
! Unchanged indices in front
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
pelem(idx:idx + da - 1) = Tens%elem(da * (ii - 1) + 1:da * ii)
end do
elseif(nc /= 0) then
! Unchanged indices in back
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
kk = ii
do jj = 1, dc
pelem(idx) = Tens%elem(kk)
idx = idx + db
kk = kk + db
end do
end do
else
! No unchanged indices in front or back
do ii = 1, db
posold(1) = posold(1) + 1
! Increment if ..
jj = 1
do while(posold(jj) >= Tens%dl(na + jj))
posold(jj) = 0
posold(jj + 1) = posold(jj + 1) + 1
jj = jj + 1
end do
! Position in the permuted array
idx = sum(posold * cumdim) + 1
pelem(idx) = Tens%elem(ii)
end do
end if
deallocate(Tens%elem)
Tens%elem => pelem
pelem => null()
!Tens%elem(:Tens%dim) = pelem
Tens%dl = Tens%dl(Tens%idx)
Tens%idx = [(ii, ii = 1, Tens%rank)]
deallocate(newdim, posold, cumdim)
!deallocate(pelem, newdim, posold, cumdim)
end subroutine permute_tensorc
"""
return
[docs]def perturb_tensor():
"""
fortran-subroutine - December 2018 (dj)
Perturb a tensor by a some small epsilon, entry-by-entry.
**Arguments**
Tens : TYPE(tensor), inout
Tensor to be perturbed.
epsilon : real, OPTIONAL, in
Scale a randomized tensor by epsilon and add to original tensor.
Default to 1e-8
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine perturb_tensor(Tens, epsilon, errst)
type(tensor), intent(inout) :: Tens
real(KIND=rKind), intent(in), optional :: epsilon
integer, intent(out), optional :: errst
! Local variables
! ---------------
! duplicate for optional arguments
real(KIND=rKind) :: delta
! random tensor
type(tensor) :: Tmp
!if(present(errst)) errst= 0
delta = 1e-8_rKind
if(present(epsilon)) delta = epsilon
call create(Tmp, Tens%dl, 'R', errst=errst)
!if(prop_error('perturb_tensor : create failed', &
! 'Tensors_include.f90:2218', errst=errst)) return
call gaxpy(Tens, delta, Tmp, errst=errst)
!if(prop_error('perturb_tensor : gaxpy failed', &
! 'Tensors_include.f90:2222', errst=errst)) return
call destroy(Tmp)
end subroutine perturb_tensor
"""
return
[docs]def perturb_tensorc():
"""
fortran-subroutine - December 2018 (dj)
Perturb a tensor by a some small epsilon, entry-by-entry.
**Arguments**
Tens : TYPE(tensorc), inout
Tensor to be perturbed.
epsilon : real, OPTIONAL, in
Scale a randomized tensor by epsilon and add to original tensor.
Default to 1e-8
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine perturb_tensorc(Tens, epsilon, errst)
type(tensorc), intent(inout) :: Tens
real(KIND=rKind), intent(in), optional :: epsilon
integer, intent(out), optional :: errst
! Local variables
! ---------------
! duplicate for optional arguments
real(KIND=rKind) :: delta
! random tensor
type(tensorc) :: Tmp
!if(present(errst)) errst= 0
delta = 1e-8_rKind
if(present(epsilon)) delta = epsilon
call create(Tmp, Tens%dl, 'R', errst=errst)
!if(prop_error('perturb_tensorc : create failed', &
! 'Tensors_include.f90:2218', errst=errst)) return
call gaxpy(Tens, delta, Tmp, errst=errst)
!if(prop_error('perturb_tensorc : gaxpy failed', &
! 'Tensors_include.f90:2222', errst=errst)) return
call destroy(Tmp)
end subroutine perturb_tensorc
"""
return
[docs]def pointto_tensor():
"""
fortran-subroutine - March 2018 (dj)
Copy Objb to Obja by setting the pointers and nullify Objb.
**Arguments**
Obja : TYPE(tensor), inout
On exit, it contains the same information as Objb
Objb : TYPE(tensor), inout
Empty on exit.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine pointto_tensor(Obja, Objb)
type(tensor), intent(inout) :: Obja, Objb
! No local variables
! ------------------
Obja%rank = Objb%rank
Objb%rank = -1
Obja%dim = Objb%dim
Objb%dim = 0
Obja%dl => Objb%dl
Objb%dl => null()
Obja%mdl => Objb%mdl
Objb%mdl => null()
Obja%idx => Objb%idx
Objb%idx => null()
Obja%elem => Objb%elem
Objb%elem => null()
end subroutine pointto_tensor
"""
return
[docs]def pointto_tensorc():
"""
fortran-subroutine - March 2018 (dj)
Copy Objb to Obja by setting the pointers and nullify Objb.
**Arguments**
Obja : TYPE(tensorc), inout
On exit, it contains the same information as Objb
Objb : TYPE(tensorc), inout
Empty on exit.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine pointto_tensorc(Obja, Objb)
type(tensorc), intent(inout) :: Obja, Objb
! No local variables
! ------------------
Obja%rank = Objb%rank
Objb%rank = -1
Obja%dim = Objb%dim
Objb%dim = 0
Obja%dl => Objb%dl
Objb%dl => null()
Obja%mdl => Objb%mdl
Objb%mdl => null()
Obja%idx => Objb%idx
Objb%idx => null()
Obja%elem => Objb%elem
Objb%elem => null()
end subroutine pointto_tensorc
"""
return
[docs]def pointto_tensorc_tensor():
"""
fortran-subroutine - March 2018 (dj)
Copy Objb to Obja by setting the pointers and nullify Objb.
Conversion from real to complex is without pointers.
**Arguments**
Obja : TYPE(tensorc), inout
On exit, it contains the same information as Objb
Objb : TYPE(tensor), inout
Empty on exit.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine pointto_tensorc_tensor(Obja, Objb)
type(tensorc), intent(inout) :: Obja
type(tensor), intent(inout) :: Objb
! No local variables
! ------------------
Obja%rank = Objb%rank
Objb%rank = -1
Obja%dim = Objb%dim
Objb%dim = 0
Obja%dl => Objb%dl
Objb%dl => null()
Obja%mdl => Objb%mdl
Objb%mdl => null()
Obja%idx => Objb%idx
Objb%idx => null()
allocate(Obja%elem(product(Obja%mdl)))
Obja%elem(:Obja%dim) = Objb%elem(:Obja%dim)
deallocate(Objb%elem)
end subroutine pointto_tensorc_tensor
"""
return
[docs]def print_tensor():
"""
fortran-subroutine - March 2016 (updated dj)
Tensor print.
**Arguments**
Obj : TYPE(tensor), in
Deprecated subroutines. Uses `write_tensor` internally.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine print_tensor(Obj)
type(tensor), intent(in) :: Obj
call write_tensor(Obj, 6, "6")
end subroutine print_tensor
"""
return
[docs]def print_tensorc():
"""
fortran-subroutine - March 2016 (updated dj)
Tensor print.
**Arguments**
Obj : TYPE(tensorc), in
Deprecated subroutines. Uses `write_tensorc` internally.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine print_tensorc(Obj)
type(tensorc), intent(in) :: Obj
call write_tensorc(Obj, 6, "6")
end subroutine print_tensorc
"""
return
[docs]def project_tensor():
"""
fortran-subroutine - March 2016 (updated dj)
Apply the projector P to Aa and store in Bb, where P is defined as
:math:`P = 1 - \sum_{\\alpha} | psiProjs_{\\alpha}> <psiProjs_{\\alpha}|`
**Arguments**
Bb : TYPE(tensor), out
Store projection of A in this tensor.
Aa : TYPE(tensor), in
Get projection of this tensor.
PsiProjs : TYPE(tensor)(\*), in
Array of tensors defining the projector.
**Details**
Used in orthogonalizing an MPS against another set of MPSs.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine project_tensor(Bb, Aa, PsiProjs, errst)
type(tensor), intent(out) :: Bb
type(tensor), intent(in) :: Aa
type(tensor), pointer, intent(in) :: PsiProjs(:)
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: jj
! overlap between tensors (dot product)
real(KIND=rKind) :: proj
!if(present(errst)) errst = 0
!do jj = 1, size(PsiProjs)
! if(Aa%rank /= PsiProjs(jj)%rank) then
! errst = raise_error('project_tensor: rank '// &
! 'mismatch.', 2, &
! 'Tensors_include.f90:2419', errst=errst)
! return
! end if
!end do
!do jj = 1, size(PsiProjs)
! if((Aa%rank /= 0) .and. any(Aa%dl /= PsiProjs(jj)%dl)) then
! errst = raise_error('project_tensor: dimension '// &
! 'mismatch.', 3, &
! 'Tensors_include.f90:2428', errst=errst)
! return
! end if
!end do
call copy(Bb, Aa)
do jj = 1, size(PsiProjs)
proj = - dot(PsiProjs(jj), Aa)
call gaxpy(Bb, proj, PsiProjs(jj))
end do
end subroutine project_tensor
"""
return
[docs]def project_tensorc():
"""
fortran-subroutine - March 2016 (updated dj)
Apply the projector P to Aa and store in Bb, where P is defined as
:math:`P = 1 - \sum_{\\alpha} | psiProjs_{\\alpha}> <psiProjs_{\\alpha}|`
**Arguments**
Bb : TYPE(tensorc), out
Store projection of A in this tensor.
Aa : TYPE(tensorc), in
Get projection of this tensor.
PsiProjs : TYPE(tensorc)(\*), in
Array of tensors defining the projector.
**Details**
Used in orthogonalizing an MPS against another set of MPSs.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine project_tensorc(Bb, Aa, PsiProjs, errst)
type(tensorc), intent(out) :: Bb
type(tensorc), intent(in) :: Aa
type(tensorc), pointer, intent(in) :: PsiProjs(:)
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: jj
! overlap between tensors (dot product)
complex(KIND=rKind) :: proj
!if(present(errst)) errst = 0
!do jj = 1, size(PsiProjs)
! if(Aa%rank /= PsiProjs(jj)%rank) then
! errst = raise_error('project_tensorc: rank '// &
! 'mismatch.', 2, &
! 'Tensors_include.f90:2419', errst=errst)
! return
! end if
!end do
!do jj = 1, size(PsiProjs)
! if((Aa%rank /= 0) .and. any(Aa%dl /= PsiProjs(jj)%dl)) then
! errst = raise_error('project_tensorc: dimension '// &
! 'mismatch.', 3, &
! 'Tensors_include.f90:2428', errst=errst)
! return
! end if
!end do
call copy(Bb, Aa)
do jj = 1, size(PsiProjs)
proj = - dot(PsiProjs(jj), Aa)
call gaxpy(Bb, proj, PsiProjs(jj))
end do
end subroutine project_tensorc
"""
return
[docs]def qmattomat_tensor():
"""
fortran-subroutine - ?? (dj)
Dummy interface resulting in a copy.
**Arguments**
Tout : TYPE(tensor), inout
Copy of input tensor.
Tin : TYPE(tensor), in
Copy this tensor, which is already a tensor without
symmetries.
Imapper : TYPE(IMAP), in
Mapping quantum numbers to matrix. Unused.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine qmattomat_tensor(Tout, Tin, Imapper, errst)
type(tensor), intent(inout) :: Tout
type(tensor), intent(in) :: Tin
type(imap), intent(in) :: Imapper
integer, intent(out), optional :: errst
! No local variables
! ------------------
call copy(Tout, Tin, errst=errst)
!if(prop_error('qmattomat_tensor : copy failed.', &
! 'Tensors_include.f90:2485', errst=errst)) return
end subroutine qmattomat_tensor
"""
return
[docs]def qmattomat_tensorc():
"""
fortran-subroutine - ?? (dj)
Dummy interface resulting in a copy.
**Arguments**
Tout : TYPE(tensorc), inout
Copy of input tensor.
Tin : TYPE(tensorc), in
Copy this tensor, which is already a tensor without
symmetries.
Imapper : TYPE(IMAP), in
Mapping quantum numbers to matrix. Unused.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine qmattomat_tensorc(Tout, Tin, Imapper, errst)
type(tensorc), intent(inout) :: Tout
type(tensorc), intent(in) :: Tin
type(imap), intent(in) :: Imapper
integer, intent(out), optional :: errst
! No local variables
! ------------------
call copy(Tout, Tin, errst=errst)
!if(prop_error('qmattomat_tensorc : copy failed.', &
! 'Tensors_include.f90:2485', errst=errst)) return
end subroutine qmattomat_tensorc
"""
return
[docs]def randomize_tensor():
"""
fortran-subroutine - Randomize the elements of the tensor.
March 2016 (updated dj)
**Arguments**
Tens : TYPE(tensor), inout
Randomize entries between [-1, 1] (real array) or [-1-i, 1+i]
(complex arrays).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine randomize_tensor(Tens, errst)
type(tensor), intent(inout) :: Tens
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
call randomize_real(Tens%elem(:Tens%dim), Tens%dim)
end subroutine randomize_tensor
"""
return
[docs]def randomize_tensorc():
"""
fortran-subroutine - Randomize the elements of the tensor.
March 2016 (updated dj)
**Arguments**
Tens : TYPE(tensorc), inout
Randomize entries between [-1, 1] (real array) or [-1-i, 1+i]
(complex arrays).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine randomize_tensorc(Tens, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
call randomize_complex(Tens%elem(:Tens%dim), Tens%dim)
end subroutine randomize_tensorc
"""
return
[docs]def rank_tensor():
"""
fortran-function - June 2017 (dj)
Return the rank of a tensor.
**Arguments**
Tens : TYPE(tensor), in
Return the rank of the tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function rank_tensor(Tens) result(myrank)
type(tensor), intent(in) :: Tens
integer :: myrank
myrank = Tens%rank
end function rank_tensor
"""
return
[docs]def rank_tensorc():
"""
fortran-function - June 2017 (dj)
Return the rank of a tensor.
**Arguments**
Tens : TYPE(tensorc), in
Return the rank of the tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function rank_tensorc(Tens) result(myrank)
type(tensorc), intent(in) :: Tens
integer :: myrank
myrank = Tens%rank
end function rank_tensorc
"""
return
[docs]def read_tensor():
"""
fortran-subroutine - August 2015 (dj)
Read an tensor from a given unit.
**Arguments**
Tens : TYPE(tensor), out
will be read from given destination.
unit : INTEGER, in
read from this unit
form : CHARACTER, in
'H' (human readable), 'B' (binary).
**Details**
For a detailed description go to the documentation of the
corresponding write method `write_tensor`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine read_tensor(Tens, unit, form, errst)
type(tensor), intent(out) :: Tens
integer, intent(in) :: unit
character, intent(in) :: form
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: i1!, i2, i3
!if(present(errst)) errst = 0
if(form == "6") then
!errst = raise_error('read_tensor: form=6 not '// &
! 'enabled!', 4, 'Tensors_include.f90:2600', errst=errst)
elseif(form == "H") then
read(unit, *) Tens%rank
allocate(Tens%dl(Tens%rank), Tens%mdl(Tens%rank), &
Tens%idx(Tens%rank))
read(unit, *) Tens%dl, Tens%mdl, Tens%idx
Tens%dim = product(Tens%dl)
allocate(Tens%elem(product(Tens%mdl)))
do i1 = 1, Tens%dim
read(unit, *) Tens%elem(i1)
end do
elseif(form =="B") then
read(unit) Tens%rank
allocate(Tens%dl(Tens%rank), Tens%mdl(Tens%rank), &
Tens%idx(Tens%rank))
read(unit) Tens%dl, Tens%mdl, Tens%idx
Tens%dim = product(Tens%dl)
allocate(Tens%elem(product(Tens%mdl)))
read(unit) Tens%elem(:Tens%dim)
else
!errst = raise_error('read_tensor: unknown '// &
! 'argument.', 4, 'Tensors_include.f90:2625', errst=errst)
end if
end subroutine read_tensor
"""
return
[docs]def read_tensorc():
"""
fortran-subroutine - August 2015 (dj)
Read an tensor from a given unit.
**Arguments**
Tens : TYPE(tensorc), out
will be read from given destination.
unit : INTEGER, in
read from this unit
form : CHARACTER, in
'H' (human readable), 'B' (binary).
**Details**
For a detailed description go to the documentation of the
corresponding write method `write_tensorc`.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine read_tensorc(Tens, unit, form, errst)
type(tensorc), intent(out) :: Tens
integer, intent(in) :: unit
character, intent(in) :: form
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: i1!, i2, i3
!if(present(errst)) errst = 0
if(form == "6") then
!errst = raise_error('read_tensorc: form=6 not '// &
! 'enabled!', 4, 'Tensors_include.f90:2600', errst=errst)
elseif(form == "H") then
read(unit, *) Tens%rank
allocate(Tens%dl(Tens%rank), Tens%mdl(Tens%rank), &
Tens%idx(Tens%rank))
read(unit, *) Tens%dl, Tens%mdl, Tens%idx
Tens%dim = product(Tens%dl)
allocate(Tens%elem(product(Tens%mdl)))
do i1 = 1, Tens%dim
read(unit, *) Tens%elem(i1)
end do
elseif(form =="B") then
read(unit) Tens%rank
allocate(Tens%dl(Tens%rank), Tens%mdl(Tens%rank), &
Tens%idx(Tens%rank))
read(unit) Tens%dl, Tens%mdl, Tens%idx
Tens%dim = product(Tens%dl)
allocate(Tens%elem(product(Tens%mdl)))
read(unit) Tens%elem(:Tens%dim)
else
!errst = raise_error('read_tensorc: unknown '// &
! 'argument.', 4, 'Tensors_include.f90:2625', errst=errst)
end if
end subroutine read_tensorc
"""
return
[docs]def scale_tensor_real():
"""
fortran-subroutine - March 2016 (updated dj)
Scale a tensor by a scalar.
**Arguments**
sc : REALSCALAR, in
Scalar to be multiplied with tensor.
Obj : TYPE(COMPLEXtensor), inout
Tensor to be scaled with scalar.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine scale_tensor_real(sc, Obj, errst)
real(KIND=rKind), intent(in) :: sc
type(tensor), intent(inout) :: Obj
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
!if(.not. associated(Obj%elem)) then
! errst = raise_error('scale_tensor_real : '//&
! 'allocation error', 99, 'Tensors_include.f90:2668', &
! errst=errst)
! return
!end if
! DSCAL / ZDSCAL / ZSCAL
call DSCAL(Obj%dim, sc, Obj%elem, 1)
end subroutine scale_tensor_real
"""
return
[docs]def scale_tensorc_complex():
"""
fortran-subroutine - March 2016 (updated dj)
Scale a tensor by a scalar.
**Arguments**
sc : REALSCALAR, in
Scalar to be multiplied with tensor.
Obj : TYPE(COMPLEXtensorc), inout
Tensor to be scaled with scalar.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine scale_tensorc_complex(sc, Obj, errst)
complex(KIND=rKind), intent(in) :: sc
type(tensorc), intent(inout) :: Obj
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
!if(.not. associated(Obj%elem)) then
! errst = raise_error('scale_tensorc_complex : '//&
! 'allocation error', 99, 'Tensors_include.f90:2668', &
! errst=errst)
! return
!end if
! DSCAL / ZDSCAL / ZSCAL
call ZSCAL(Obj%dim, sc, Obj%elem, 1)
end subroutine scale_tensorc_complex
"""
return
[docs]def scale_tensorc_real():
"""
fortran-subroutine - March 2016 (updated dj)
Scale a tensor by a scalar.
**Arguments**
sc : REALSCALAR, in
Scalar to be multiplied with tensor.
Obj : TYPE(COMPLEXtensorc), inout
Tensor to be scaled with scalar.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine scale_tensorc_real(sc, Obj, errst)
real(KIND=rKind), intent(in) :: sc
type(tensorc), intent(inout) :: Obj
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
!if(.not. associated(Obj%elem)) then
! errst = raise_error('scale_tensorc_real : '//&
! 'allocation error', 99, 'Tensors_include.f90:2668', &
! errst=errst)
! return
!end if
! DSCAL / ZDSCAL / ZSCAL
call ZDSCAL(Obj%dim, sc, Obj%elem, 1)
end subroutine scale_tensorc_real
"""
return
[docs]def set_hash_tensor():
"""
fortran-subroutine - October 2016 (dj)
Dummy subroutine to comply with interface in qTensors. No effect
on usual tensor.
**Arguments**
Tens : TYPE(tensor), inout
Tens is returned without modification as it has no symmetries.
idxs : INTEGER(\*), in
Specify indices to be hashed in regard to the legs of the tensor.
ii : INTEGER, OPTIONAL, in
Specify index in terms of array of tensors for generating new hash.
If not present, all tensors are set again.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine set_hash_tensor(Tens, idxs, ii, errst)
type(tensor), intent(inout) :: Tens
integer, intent(in) :: idxs(:)
integer, intent(in), optional :: ii
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
end subroutine set_hash_tensor
"""
return
[docs]def set_hash_tensorc():
"""
fortran-subroutine - October 2016 (dj)
Dummy subroutine to comply with interface in qTensors. No effect
on usual tensor.
**Arguments**
Tens : TYPE(tensorc), inout
Tens is returned without modification as it has no symmetries.
idxs : INTEGER(\*), in
Specify indices to be hashed in regard to the legs of the tensor.
ii : INTEGER, OPTIONAL, in
Specify index in terms of array of tensors for generating new hash.
If not present, all tensors are set again.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine set_hash_tensorc(Tens, idxs, ii, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: idxs(:)
integer, intent(in), optional :: ii
integer, intent(out), optional :: errst
!if(present(errst)) errst = 0
end subroutine set_hash_tensorc
"""
return
[docs]def split_tensor():
"""
fortran-subroutine - October 2016 (dj)
Split one index of a tensor into multiple indices.
**Arguments**
Tens : TYPE(tensor), inout
Split one index of this tensor.
idx : INTEGER, in
Split this index into multiple indices.
Sl : TYPE(splitlink), in
Information about splitting the tensor.
sidx : INTEGER, OPTIONAL, in
Specify which information in Sl is accessed. Default to 1.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine split_tensor(Tens, idx, Sl, sidx, errst)
type(tensor), intent(inout) :: Tens
integer, intent(in) :: idx
type(splitlink), intent(in) :: Sl
integer, intent(in), optional :: sidx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! old rank
integer :: orank
! index considering permutation
integer :: kk
! copy of dimensions / permutations
integer, dimension(:), pointer :: darr
! new dimensions of the splitted link
integer, dimension(:), allocatable :: dims
!if(present(errst)) errst = 0
! Copy the dimensions
allocate(dims(size(Sl%dims, 1)))
if(present(sidx)) then
dims = Sl%dims(:, sidx)
else
dims = Sl%dims(:, 1)
end if
call permute(Tens, errst=errst)
!if(prop_error('split_tensor: permute failed.', &
! 'Tensors_include.f90:2795', errst=errst)) return
kk = idx
!kk = Tens%idx(idx)
!if(product(dims) /= Tens%dl(kk)) then
! errst = raise_error('split_tensor: dimension '// &
! 'mismatch.', 3, 'Tensors_include.f90:2802', errst=errst)
! return
!end if
allocate(darr(Tens%rank))
darr = Tens%dl
orank = Tens%rank
Tens%rank = Tens%rank + size(dims) - 1
deallocate(Tens%dl, Tens%idx)
allocate(Tens%dl(Tens%rank), Tens%idx(Tens%rank))
Tens%idx = [(ii, ii = 1, Tens%rank)]
! Reset dimensions
Tens%dl(:kk - 1) = darr(:kk - 1)
Tens%dl(kk:kk + size(dims) - 1) = dims
Tens%dl(kk + size(dims):) = darr(kk + 1:)
deallocate(darr, dims)
end subroutine split_tensor
"""
return
[docs]def split_all_tensor():
"""
fortran-subroutine - October 2017 (dj)
Split all links of a tensor.
**Arguments**
Tens : TYPE(tensor), inout
Split all links of this tensor.
Sls : TYPE(splitlink)(\*), in
Contains the information how to split each index. Dimension
of Sls must correspond to the rank of the tensor.
sidx : INTEGER, OPTIONAL, in
Choose an index within the splitlink information. Default to 1.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine split_all_tensor(Tens, Sls, sidx, errst)
type(tensor), intent(inout) :: Tens
type(splitlink), dimension(:), intent(in) :: Sls
integer, dimension(:), intent(in), optional :: sidx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! duplette for sidx variable
integer :: kk
!if(present(errst)) errst = 0
do ii = size(Sls), 1, -1
kk = 1
if(present(sidx)) kk = sidx(ii)
call split(Tens, ii, Sls(ii), sidx=kk, errst=errst)
!if(prop_error('split_all_tensor : split failed.', &
! 'Tensors_include.f90:2876', errst=errst)) return
end do
end subroutine split_all_tensor
"""
return
[docs]def split_tensorc():
"""
fortran-subroutine - October 2016 (dj)
Split one index of a tensor into multiple indices.
**Arguments**
Tens : TYPE(tensorc), inout
Split one index of this tensor.
idx : INTEGER, in
Split this index into multiple indices.
Sl : TYPE(splitlink), in
Information about splitting the tensor.
sidx : INTEGER, OPTIONAL, in
Specify which information in Sl is accessed. Default to 1.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine split_tensorc(Tens, idx, Sl, sidx, errst)
type(tensorc), intent(inout) :: Tens
integer, intent(in) :: idx
type(splitlink), intent(in) :: Sl
integer, intent(in), optional :: sidx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! old rank
integer :: orank
! index considering permutation
integer :: kk
! copy of dimensions / permutations
integer, dimension(:), pointer :: darr
! new dimensions of the splitted link
integer, dimension(:), allocatable :: dims
!if(present(errst)) errst = 0
! Copy the dimensions
allocate(dims(size(Sl%dims, 1)))
if(present(sidx)) then
dims = Sl%dims(:, sidx)
else
dims = Sl%dims(:, 1)
end if
call permute(Tens, errst=errst)
!if(prop_error('split_tensorc: permute failed.', &
! 'Tensors_include.f90:2795', errst=errst)) return
kk = idx
!kk = Tens%idx(idx)
!if(product(dims) /= Tens%dl(kk)) then
! errst = raise_error('split_tensorc: dimension '// &
! 'mismatch.', 3, 'Tensors_include.f90:2802', errst=errst)
! return
!end if
allocate(darr(Tens%rank))
darr = Tens%dl
orank = Tens%rank
Tens%rank = Tens%rank + size(dims) - 1
deallocate(Tens%dl, Tens%idx)
allocate(Tens%dl(Tens%rank), Tens%idx(Tens%rank))
Tens%idx = [(ii, ii = 1, Tens%rank)]
! Reset dimensions
Tens%dl(:kk - 1) = darr(:kk - 1)
Tens%dl(kk:kk + size(dims) - 1) = dims
Tens%dl(kk + size(dims):) = darr(kk + 1:)
deallocate(darr, dims)
end subroutine split_tensorc
"""
return
[docs]def split_all_tensorc():
"""
fortran-subroutine - October 2017 (dj)
Split all links of a tensor.
**Arguments**
Tens : TYPE(tensorc), inout
Split all links of this tensor.
Sls : TYPE(splitlink)(\*), in
Contains the information how to split each index. Dimension
of Sls must correspond to the rank of the tensor.
sidx : INTEGER, OPTIONAL, in
Choose an index within the splitlink information. Default to 1.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine split_all_tensorc(Tens, Sls, sidx, errst)
type(tensorc), intent(inout) :: Tens
type(splitlink), dimension(:), intent(in) :: Sls
integer, dimension(:), intent(in), optional :: sidx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
! duplette for sidx variable
integer :: kk
!if(present(errst)) errst = 0
do ii = size(Sls), 1, -1
kk = 1
if(present(sidx)) kk = sidx(ii)
call split(Tens, ii, Sls(ii), sidx=kk, errst=errst)
!if(prop_error('split_all_tensorc : split failed.', &
! 'Tensors_include.f90:2876', errst=errst)) return
end do
end subroutine split_all_tensorc
"""
return
[docs]def trace_tensor():
"""
fortran-function - June 2017 (dj)
Returns the trace for rank-2 tensors.
**Arguments**
Tens : TYPE(tensor), in
Calculate the trace assuming it is a rank-2 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function trace_tensor(Tens, errst) result(tr)
type(tensor), intent(in) :: Tens
integer, intent(out), optional :: errst
real(KIND=rKind) :: tr
! External functions
! ------------------
! goes into DDOT,ZDOTC
!real(KIND=rKind), external :: dsum
! Local variables
! ---------------
! for looping
integer :: ii
! next index
integer :: idx
!CHEK_if(present(errst)) errst = 0
!if(Tens%rank /= 2) then
! errst = raise_error('trace_tensor: not rank 2.', &
! 99, 'Tensors_include.f90:2926', errst=errst)
! return
!end if
!tr = dsum(min(Tens%dl(1), Tens%dl(2)), Tens%elem, Tens%dl(1) + 1)
idx = 1
tr = Tens%elem(idx)
do ii = 2, Tens%dl(1)
idx = idx + Tens%dl(1) + 1
tr = tr + Tens%elem(idx)
end do
end function trace_tensor
"""
return
[docs]def trace_tensorc():
"""
fortran-function - June 2017 (dj)
Returns the trace for rank-2 tensors.
**Arguments**
Tens : TYPE(tensorc), in
Calculate the trace assuming it is a rank-2 tensor.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
function trace_tensorc(Tens, errst) result(tr)
type(tensorc), intent(in) :: Tens
integer, intent(out), optional :: errst
complex(KIND=rKind) :: tr
! External functions
! ------------------
! goes into DDOT,ZDOTC
!complex(KIND=rKind), external :: zsum
! Local variables
! ---------------
! for looping
integer :: ii
! next index
integer :: idx
!CHEK_if(present(errst)) errst = 0
!if(Tens%rank /= 2) then
! errst = raise_error('trace_tensorc: not rank 2.', &
! 99, 'Tensors_include.f90:2926', errst=errst)
! return
!end if
!tr = zsum(min(Tens%dl(1), Tens%dl(2)), Tens%elem, Tens%dl(1) + 1)
idx = 1
tr = Tens%elem(idx)
do ii = 2, Tens%dl(1)
idx = idx + Tens%dl(1) + 1
tr = tr + Tens%elem(idx)
end do
end function trace_tensorc
"""
return
[docs]def transposed_tensor():
"""
fortran-subroutine - July 2016 (dj)
Transposition of indices or permutation.
**Arguments**
Tens : TYPE(tensor), inout
Save a transposition/permutation on the indices of this tensor.
perm : INTEGER(\*), OPTIONAL, in
permutation array has length equal to the rank of the tensor with
unique entries 1 to rank.
Default to rank, rank - 1, ..., 2, 1 (transpose)
doperm : LOGICAL, OPTIONAL, in
If true, then permute immediately. If false or not present,
permutation is only indicated in the indices.
**Details**
The permutation is only stored in the indices and no actions
on the actual memory are carried out. The permutation array [3, 1, 2]
is to be read as "The new index 1 comes from old index 3" etc.
(template defined in Tensors_include.f90)
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine transposed_tensor(Tens, perm, doperm, errst)
type(tensor), intent(inout) :: Tens
integer, dimension(:), intent(in), optional :: perm
logical, intent(in), optional :: doperm
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
!if(present(errst)) errst = 0
! Fast return
if(Tens%rank < 2) return
! Default transpose - reverse indices
if(.not. present(perm)) then
Tens%idx = [(ii, ii = Tens%rank, 1, -1)]
else
!if(size(perm) /= Tens%rank) then
! errst = raise_error('transpose_tensor: rank '// &
! 'mismatch.', 2, 'Tensors_include.f90:3007', errst=errst)
! return
!end if
Tens%idx = Tens%idx(perm)
end if
if(present(doperm)) then
if(doperm) then
call permute(Tens, errst=errst)
!if(prop_error('transposed_tensor: permute '//&
! 'failed.', 'Tensors_include.f90:3018', errst=errst)) return
end if
end if
end subroutine transposed_tensor
"""
return
[docs]def write_tensor():
"""
fortran-subroutine - August 2015 (dj)
Write the tensor to a file or the std-out.
**Arguments**
Tens : TYPE(tensor), in
will be written to given destination.
unit : INTEGER, in
write on this unit
form : CHARACTER, in
'H' (human readable), 'B' (binary) or '6' (non-zero values
intended for standard-output).
**Details**
For form = `H` or `B` the information about the rank-3 tensor is
written in the following form:
1) three integers, dimension of the tensor
2) `H` : Entries of the matrix where each entry is in its own line. The
tensor is saved row-major (meaning the order of the elements is
(1,1,1), (1,1,2), (1,1,3), ..., (1,1,d3), (1,2,1), (1,2,2), ...).
`B` : elements are written in one step (column-major!).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine write_tensor(Tens, unit, form, errst)
TYPE(tensor), intent(in) :: Tens
integer, intent(in) :: unit
character, intent(in) :: form
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: i1!, i2, i3
!if(present(errst)) errst = 0
if(form == "6") then
write(unit, *) 'Tensor rank', Tens%rank
write(unit, *) 'Tensor dimensions', Tens%dl, Tens%mdl, Tens%idx
do i1 = 1, Tens%dim
if(abs(Tens%elem(i1)) > 10.0_rKind**(-6)) then
write(unit, *) 'i, M(i)', i1, Tens%elem(i1)
end if
end do
elseif(form == "H") then
write(unit, *) Tens%rank
write(unit, *) Tens%dl, Tens%mdl, Tens%idx
do i1 = 1, Tens%dim
write(unit, *) Tens%elem(i1)
end do
elseif(form =="B") then
write(unit) Tens%rank
write(unit) Tens%dl, Tens%mdl, Tens%idx
write(unit) Tens%elem(:Tens%dim)
else
!errst = raise_error('write_tensor: unknown argument.', &
! 4, 'Tensors_include.f90:3100', errst=errst)
end if
end subroutine write_tensor
"""
return
[docs]def transposed_tensorc():
"""
fortran-subroutine - July 2016 (dj)
Transposition of indices or permutation.
**Arguments**
Tens : TYPE(tensorc), inout
Save a transposition/permutation on the indices of this tensor.
perm : INTEGER(\*), OPTIONAL, in
permutation array has length equal to the rank of the tensor with
unique entries 1 to rank.
Default to rank, rank - 1, ..., 2, 1 (transpose)
doperm : LOGICAL, OPTIONAL, in
If true, then permute immediately. If false or not present,
permutation is only indicated in the indices.
**Details**
The permutation is only stored in the indices and no actions
on the actual memory are carried out. The permutation array [3, 1, 2]
is to be read as "The new index 1 comes from old index 3" etc.
(template defined in Tensors_include.f90)
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine transposed_tensorc(Tens, perm, doperm, errst)
type(tensorc), intent(inout) :: Tens
integer, dimension(:), intent(in), optional :: perm
logical, intent(in), optional :: doperm
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: ii
!if(present(errst)) errst = 0
! Fast return
if(Tens%rank < 2) return
! Default transpose - reverse indices
if(.not. present(perm)) then
Tens%idx = [(ii, ii = Tens%rank, 1, -1)]
else
!if(size(perm) /= Tens%rank) then
! errst = raise_error('transpose_tensorc: rank '// &
! 'mismatch.', 2, 'Tensors_include.f90:3007', errst=errst)
! return
!end if
Tens%idx = Tens%idx(perm)
end if
if(present(doperm)) then
if(doperm) then
call permute(Tens, errst=errst)
!if(prop_error('transposed_tensorc: permute '//&
! 'failed.', 'Tensors_include.f90:3018', errst=errst)) return
end if
end if
end subroutine transposed_tensorc
"""
return
[docs]def write_tensorc():
"""
fortran-subroutine - August 2015 (dj)
Write the tensor to a file or the std-out.
**Arguments**
Tens : TYPE(tensorc), in
will be written to given destination.
unit : INTEGER, in
write on this unit
form : CHARACTER, in
'H' (human readable), 'B' (binary) or '6' (non-zero values
intended for standard-output).
**Details**
For form = `H` or `B` the information about the rank-3 tensor is
written in the following form:
1) three integers, dimension of the tensor
2) `H` : Entries of the matrix where each entry is in its own line. The
tensor is saved row-major (meaning the order of the elements is
(1,1,1), (1,1,2), (1,1,3), ..., (1,1,d3), (1,2,1), (1,2,2), ...).
`B` : elements are written in one step (column-major!).
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine write_tensorc(Tens, unit, form, errst)
TYPE(tensorc), intent(in) :: Tens
integer, intent(in) :: unit
character, intent(in) :: form
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for looping
integer :: i1!, i2, i3
!if(present(errst)) errst = 0
if(form == "6") then
write(unit, *) 'Tensor rank', Tens%rank
write(unit, *) 'Tensor dimensions', Tens%dl, Tens%mdl, Tens%idx
do i1 = 1, Tens%dim
if(abs(Tens%elem(i1)) > 10.0_rKind**(-6)) then
write(unit, *) 'i, M(i)', i1, Tens%elem(i1)
end if
end do
elseif(form == "H") then
write(unit, *) Tens%rank
write(unit, *) Tens%dl, Tens%mdl, Tens%idx
do i1 = 1, Tens%dim
write(unit, *) Tens%elem(i1)
end do
elseif(form =="B") then
write(unit) Tens%rank
write(unit) Tens%dl, Tens%mdl, Tens%idx
write(unit) Tens%elem(:Tens%dim)
else
!errst = raise_error('write_tensorc: unknown argument.', &
! 4, 'Tensors_include.f90:3100', errst=errst)
end if
end subroutine write_tensorc
"""
return
[docs]def write_as_matrix_tensor():
"""
fortran-subroutine - August 2017 (dj, updated)
Write the tensor as a matrix in column-major format. If complex,
real part of all entries before complex part of all entries
**Arguments**
Tens : TYPE(tensor), in
will be written to given destination.
unit : INTEGER, in
write on this unit
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine write_as_matrix_tensor(Mat, unit, errst)
type(tensor), intent(in) :: Mat
integer, intent(in) :: unit
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for building the formatting string
character(16) :: iwstring, specstring
!if(present(errst)) errst = 0
!if(Mat%rank /= 2) then
! errst = raise_error('write_as_matrix_tensor: '//&
! 'rank mismatch.', 99, &
! 'Tensors_include.f90:3148', errst=errst)
! return
!end if
write(iwstring, '(I4)') Mat%dim
specstring = "("//trim(adjustl(iwstring))//"E30.15E3)"
write(unit, specstring) real(Mat%elem(:Mat%dim), KIND=rKind)
!write(unit, specstring) aimag(Mat%elem(:Mat%dim))
end subroutine write_as_matrix_tensor
"""
return
[docs]def write_as_matrix_tensorc():
"""
fortran-subroutine - August 2017 (dj, updated)
Write the tensor as a matrix in column-major format. If complex,
real part of all entries before complex part of all entries
**Arguments**
Tens : TYPE(tensorc), in
will be written to given destination.
unit : INTEGER, in
write on this unit
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine write_as_matrix_tensorc(Mat, unit, errst)
type(tensorc), intent(in) :: Mat
integer, intent(in) :: unit
integer, intent(out), optional :: errst
! Local variables
! ---------------
! for building the formatting string
character(16) :: iwstring, specstring
!if(present(errst)) errst = 0
!if(Mat%rank /= 2) then
! errst = raise_error('write_as_matrix_tensorc: '//&
! 'rank mismatch.', 99, &
! 'Tensors_include.f90:3148', errst=errst)
! return
!end if
write(iwstring, '(I4)') Mat%dim
specstring = "("//trim(adjustl(iwstring))//"E30.15E3)"
write(unit, specstring) real(Mat%elem(:Mat%dim), KIND=rKind)
write(unit, specstring) aimag(Mat%elem(:Mat%dim))
end subroutine write_as_matrix_tensorc
"""
return
[docs]def copy_splitlink():
"""
fortran-subroutine - October 2017 (dj)
Copy splitlink object.
**Arguments**
Slout : TYPE(splitlink), inout
Copy the entries of `Slin` here.
Slin : TYPE(splitlink), in
source for copy into `Slout`
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine copy_splitlink(Slout, Slin, errst)
type(splitlink), intent(inout) :: Slout
type(splitlink), intent(in) :: Slin
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
allocate(Slout%hashes(size(Slin%hashes, 1)), &
Slout%dims(size(Slin%dims, 1), size(Slin%dims, 2)))
Slout%hashes = Slin%hashes
Slout%dims = Slin%dims
end subroutine copy_splitlink
"""
return
[docs]def create_splitlink_dimensions():
"""
fortran-subroutine - October 2017 (dj)
Create a simple splitlink passing the dimensions.
**Arguments**
Sl : TYPE(splitlink), inout
Setup this splitlink object.
dims : INTEGER(\*), in
Dimension for splitting a tensor link.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_splitlink_dimensions(Sl, dims)
type(splitlink), intent(inout) :: Sl
integer, dimension(:), intent(in) :: dims
! No local variables
! ------------------
allocate(Sl%hashes(1), Sl%dims(size(dims, 1), 1))
Sl%hashes = 0.0_rKind
Sl%dims(:, 1) = dims
end subroutine create_splitlink_dimensions
"""
return
[docs]def create_splitlink_tensor():
"""
fortran-subroutine - October 2017 (dj)
Create a splitlink object from the identity matrix.
**Arguments**
Sl : TYPE(splitlink), inout
Fill this splitlink object with the splitting for an
identity matrix.
Idop : TYPE(tensor), in
The dimension is taken from the identity to create
the splitlink object.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_splitlink_tensor(Sl, Idop)
type(splitlink), intent(inout) :: Sl
type(tensor), intent(in) :: Idop
! No local variables
! ------------------
allocate(Sl%hashes(1), Sl%dims(2, 1))
Sl%hashes = 0.0_rKind
Sl%dims(:, 1) = Idop%dl
end subroutine create_splitlink_tensor
"""
return
[docs]def create_splitlink_tensorc():
"""
fortran-subroutine - October 2017 (dj)
Create a splitlink object from the identity matrix.
**Arguments**
Sl : TYPE(splitlink), inout
Fill this splitlink object with the splitting for an
identity matrix.
Idop : TYPE(tensorc), in
The dimension is taken from the identity to create
the splitlink object.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine create_splitlink_tensorc(Sl, Idop)
type(splitlink), intent(inout) :: Sl
type(tensorc), intent(in) :: Idop
! No local variables
! ------------------
allocate(Sl%hashes(1), Sl%dims(2, 1))
Sl%hashes = 0.0_rKind
Sl%dims(:, 1) = Idop%dl
end subroutine create_splitlink_tensorc
"""
return
[docs]def destroy_splitlink():
"""
fortran-subroutine - September 2017 (dj)
Destroy the object taking care of splitting a link.
**Arguments**
Sl : TYPE(splitlink), inout
Deallocate on exit.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine destroy_splitlink(Sl)
type(splitlink), intent(inout) :: Sl
! No local variables
! ------------------
deallocate(Sl%hashes, Sl%dims)
end subroutine destroy_splitlink
"""
return
[docs]def init_splitlink_tensor():
"""
fortran-subroutine - September 2017 (dj)
Store the information about fusing links in a tensor.
**Arguments**
Sl : TYPE(splitlink), inout
Splitlink information is initialized after this subroutine.
It needs to be finalized.
Tens : TYPE(tensor), inout
The tensor to be fused.
idx : INTEGER(\*), inout
The links to be fused.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine init_splitlink_tensor(Sl, Tens, idx, errst)
type(splitlink), intent(inout) :: Sl
type(tensor), intent(in) :: Tens
integer, dimension(:), intent(in) :: idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! number of links fused
integer :: nn
!if(present(errst)) errst = 0
nn = size(idx, 1)
allocate(Sl%hashes(1), Sl%dims(nn, 1))
Sl%hashes = 0.0_rKind
Sl%dims(:, 1) = Tens%dl(idx)
end subroutine init_splitlink_tensor
"""
return
[docs]def init_splitlink_tensorc():
"""
fortran-subroutine - September 2017 (dj)
Store the information about fusing links in a tensor.
**Arguments**
Sl : TYPE(splitlink), inout
Splitlink information is initialized after this subroutine.
It needs to be finalized.
Tens : TYPE(tensorc), inout
The tensor to be fused.
idx : INTEGER(\*), inout
The links to be fused.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine init_splitlink_tensorc(Sl, Tens, idx, errst)
type(splitlink), intent(inout) :: Sl
type(tensorc), intent(in) :: Tens
integer, dimension(:), intent(in) :: idx
integer, intent(out), optional :: errst
! Local variables
! ---------------
! number of links fused
integer :: nn
!if(present(errst)) errst = 0
nn = size(idx, 1)
allocate(Sl%hashes(1), Sl%dims(nn, 1))
Sl%hashes = 0.0_rKind
Sl%dims(:, 1) = Tens%dl(idx)
end subroutine init_splitlink_tensorc
"""
return
[docs]def finalize_splitlink_tensor():
"""
fortran-subroutine - October 2017 (dj)
Finalize splitlink is a dummy function for tensors without symmetry.
**Arguments**
Sl : TYPE(splitlink), inout
On exit, complete splitlink object. On entry, the hashes are
filled with the indices of the corresponding block.
Tens : TYPE(tensor), in
The tensor which was fused.
idx : INTEGER, in
Index of the new link, which has been fused.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine finalize_splitlink_tensor(Sl, Tens, idx, errst)
type(splitlink), intent(inout) :: Sl
type(tensor), intent(in) :: Tens
integer, intent(in) :: idx
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
! Avoid unused dummy arguments for intentionally unused arguments
!if(.false.) print *, Sl%hashes
!if(.false.) print *, Tens%rank
!if(.false.) print *, idx
end subroutine finalize_splitlink_tensor
"""
return
[docs]def finalize_splitlink_tensorc():
"""
fortran-subroutine - October 2017 (dj)
Finalize splitlink is a dummy function for tensors without symmetry.
**Arguments**
Sl : TYPE(splitlink), inout
On exit, complete splitlink object. On entry, the hashes are
filled with the indices of the corresponding block.
Tens : TYPE(tensorc), in
The tensor which was fused.
idx : INTEGER, in
Index of the new link, which has been fused.
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine finalize_splitlink_tensorc(Sl, Tens, idx, errst)
type(splitlink), intent(inout) :: Sl
type(tensorc), intent(in) :: Tens
integer, intent(in) :: idx
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
! Avoid unused dummy arguments for intentionally unused arguments
!if(.false.) print *, Sl%hashes
!if(.false.) print *, Tens%rank
!if(.false.) print *, idx
end subroutine finalize_splitlink_tensorc
"""
return
[docs]def read_splitlink():
"""
fortran-subroutine - November 2017 (dj)
Write a splitlink object to an open file specified via the unit.
**Arguments**
Sl : TYPE(splitlink), in
On exit, splitlink object is filled with data from file.
unit : INTEGER, in
Unit where the file is open.
form : CHARACTER, in
Specifies type of read between human readable output ('H')
and binary in an unformatted file ('B').
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine read_splitlink(Sl, unit, form, errst)
type(splitlink), intent(inout) :: Sl
integer, intent(in) :: unit
character, intent(in) :: form
integer, intent(out), optional :: errst
! Local variables
! ---------------
! dimensions
integer :: d1, d2
if(form == 'H') then
! Human readable file
read(unit, *) d1, d2
allocate(Sl%hashes(d2), Sl%dims(d1, d2))
read(unit, *) Sl%hashes
read(unit, *) Sl%dims
elseif(form == 'B') then
! Binary file
read(unit) d1, d2
allocate(Sl%hashes(d2), Sl%dims(d1, d2))
read(unit) Sl%hashes
read(unit) Sl%dims
else
!errst = raise_error('read_splitlink : wrong formatting.', &
! 99, 'Tensors_include.f90:3481', errst=errst)
return
end if
end subroutine read_splitlink
"""
return
[docs]def write_splitlink():
"""
fortran-subroutine - November 2017 (dj)
Write a splitlink object to an open file specified via the unit.
**Arguments**
Sl : TYPE(splitlink), in
Splitlink object to be written.
unit : INTEGER, in
Unit where the file is open.
form : CHARACTER, in
Specifies type of print between debugging output ('6'), human
readable output ('H') and binary in an unformatted file ('B').
**Source Code**
.. hidden-code-block:: fortran
:label: show / hide f90 code
subroutine write_splitlink(Sl, unit, form, errst)
type(splitlink), intent(in) :: Sl
integer, intent(in) :: unit
character, intent(in) :: form
integer, intent(out), optional :: errst
! No local variables
! ------------------
!if(present(errst)) errst = 0
if(form == '6') then
write(unit, *) 'Splitlink object - hashes and dimensions:'
write(unit, *) Sl%hashes
write(unit, *) Sl%dims
elseif(form == 'H') then
! Human readable output
write(unit, *) size(Sl%dims, 1), size(Sl%dims, 2)
write(unit, *) Sl%hashes
write(unit, *) Sl%dims
elseif(form == 'B') then
! binary file
write(unit) size(Sl%dims, 1), size(Sl%dims, 2)
write(unit) Sl%hashes
write(unit) Sl%dims
else
!errst = raise_error('write_splitlink : wrong formatting.', &
! 99, 'Tensors_include.f90:3551', errst=errst)
return
end if
end subroutine write_splitlink
"""
return