add option to disable single-precision operators

fixes problems for LS-DYNA users as reported in #10. Now the include statement is done via a preprocessor flag `#include "ttb/ttb_library.F"` and all *_r4 procedures may be excluded from compiling with `#define NOR4`. Updated docs too.

README.md

@@ -35,7 +35,7 @@ The idea is to create derived data types for rank 1, rank 2 and rank 4 tensors (
 The most basic example on how to use this module is to [download the module](https://github.com/adtzlr/ttb/archive/master.zip), put the 'ttb'-Folder in your working directory and add two lines of code:
 
 ```fortran
-       include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
        program script101_ttb
        use Tensor

docs/example.md

@@ -29,7 +29,7 @@ with the fourth order identity tensor
 The two equations are now implemented in a Total Lagrange user subroutine  with the help of this Tensor module as follows:
 
 ```fortran
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/example_neohooke.md

@@ -29,7 +29,7 @@ with the fourth order identity tensor
 The two equations are now implemented in a Total Lagrange user subroutine  with the help of this Tensor module as follows:
 
 ```fortran
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/example_stvenantkirchhoff.md

@@ -15,7 +15,7 @@ and
 Before we are able to add our own user code, we have to start with an empty fortran subroutine header for MSC.Marc's HYPELA2. Similar headers are provided for Abaqus, ANSYS, etc in the corresponding manuals.
 
 ```fortran
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/examples/hypela2_nh_ttb.f

@@ -1,4 +1,4 @@
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/examples/hypela2_nh_ttb_simple.f

@@ -1,4 +1,4 @@
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/examples/hypela2_nonlinear_viscoelasticity.f

@@ -1,4 +1,4 @@
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/examples/hypela2_stvenantkirchhoff.f

@@ -1,4 +1,4 @@
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       subroutine hypela2(d,g,e,de,s,t,dt,ngens,m,nn,kcus,matus,ndi,
      2             nshear,disp,dispt,coord,ffn,frotn,strechn,eigvn,ffn1,

docs/examples/script_umat.f

@@ -1,4 +1,4 @@
-       include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
        program script_tensortoolbox
 

docs/examples/umat_nh_ttb_simple.f

@@ -1,4 +1,4 @@
-      include 'ttb/ttb_library.f'
+#include "ttb/ttb_library.F"
 
       SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,
      1 RPL,DDSDDT,DRPLDE,DRPLDT,

docs/index.md

@@ -20,12 +20,12 @@ Commercial FEM software packages offer interfaces (user subroutines written in F
 
 ## Author
 Andreas Dutzler, Graz, Austria.
-I'm a PhD student at the Insitute of Machine Components and Methods of Development, Area of Structural Durability and Railway Engineering at Graz University of Technology.
+I'm a PhD student at the Insitute of Structural Durability and Railway Technology at Graz University of Technology.
 
 ## Discussion
 For questions I opened a [Thread at PolymerFEM.com](https://polymerfem.com/forum/polymerfem-downloads/user-material-subroutines/25668-discussion-tensor-toolbox-for-modern-fortran-ttb) for discussions regarding this module.
 
 ## Credits
 Naumann, C.: [Chemisch-mechanisch gekoppelte Modellierung und Simulation oxidativer Alterungsvorgänge in Gummibauteilen (German)](http://nbn-resolving.de/urn:nbn:de:bsz:ch1-qucosa-222075). PhD thesis. Fakultät für Maschinenbau der Technischen Universität Chemnitz, 2016.
 
-## Last updated: 2018/04/06
+## Last updated: 2020/06/02

docs/installation.md

@@ -1,7 +1,15 @@
 This Toolbox is a fortran module which can be used inside modern Fortran compilers. For compatibility reasons to commercial FEM Software packages it is written in _Fixed Format_.
 
 ## Installation Requirements
-This Toolbox is tested on Windows with both Intel Fortran 2015 and GFortran 6.3. If you are using Linux it **should** work (but it is untested).
+This Toolbox is tested on Windows with both Intel Fortran >2015 (in combination with MSC.Marc) and GFortran >6.3. If you are using Linux it **should** work (but it is untested).
 
 ## Download
 [Download the module](https://github.com/adtzlr/ttb/archive/master.zip), put the `ttb`-Folder in your working directory and you are ready to dive into comfortable tensor manipulations in Fortran.
+
+## A note on LS-DYNA Users
+If you have problems as reported [here] (https://github.com/adtzlr/ttb/issues/10), please add the following line **before** the Tensor-Toolbox include statement. This deactivates tensor with single-precision scalar multiplications and divisions. **Warning**: Now take care to only use double-precision constants in your code!
+
+```fortran
+#define NOR4
+#include "ttb/ttb_library.F"
+```

docs/quickstartguide.md

@@ -2,7 +2,7 @@
 The most basic example on how to use this module is to [download the module](https://github.com/adtzlr/ttb/archive/master.zip), put the 'ttb'-Folder in your working directory and add two lines of code:
 
 ```fortran
-       include 'ttb/ttb_library.f'
+#include 'ttb/ttb_library.F'
 
        program script101_ttb
        use Tensor
@@ -12,6 +12,6 @@ The most basic example on how to use this module is to [download the module](htt
 
        end program script101_ttb
 ```
-The `include 'ttb/ttb_library.f'` statement replaces the line with the content of the ttb-module. The first line in a program or subroutine is now a `use Tensor` statement. That's it - now you're ready to go.
+The `#include 'ttb/ttb_library.F'` statement replaces the line with the content of the ttb-module. The first line in a program or subroutine is now a `use Tensor` statement. That's it - now you're ready to go. Be sure to save your files with **UPPERCASE** file endings, e.g. `file.F` instead of `file.f`. This tells the Fortran compiler to use a preprocessor.
 
 Continue to [Example](example.md) section. For a list and detailed information of available functions go [here](functions.md).

ttb/libderivativeexp.f

@@ -26,7 +26,7 @@ function dexp_2(T)
         dexp_2 = Eye.cdya.Eye
 
         I1 = T**Eye
-        I2 = 0.5*(I1**2-T**T)
+        I2 = 5.0d-1*(I1**2-T**T)
         I3 = det(T)
 
         N_i = Eye
@@ -81,7 +81,7 @@ function dexp_2s(T)
         dexp_2s = Eye.cdya.Eye
 
         I1 = T**Eye
-        I2 = 0.5*(I1**2-T**T)
+        I2 = 5.0d-1*(I1**2-T**T)
         I3 = det(T)
 
         N_i = Eye

ttb/libexp.f

@@ -22,10 +22,10 @@ function exp_2(T)
         i = 3
 
         Eye = identity2(Eye)
-        exp_2 = Eye + T + 0.5 * (T*T)
+        exp_2 = Eye + T + 5.0d-1 * (T*T)
 
         I1 = T**Eye
-        I2 = 0.5*(I1**2-T**T)
+        I2 = 5.0d-1*(I1**2-T**T)
         I3 = det(T)
 
         N_i = Eye

ttb/libinv.f

@@ -3,19 +3,19 @@ function inv_2(T)
 
         type(Tensor2), intent(in) :: T
         type(Tensor2) :: inv_2
-        real(kind=8) :: detT
+        real(kind=8) :: idetT
 
-        detT = det(T)
+        idetT = 1.d0/det(T)
 
-        inv_2%ab(1,1)=+1./detT*(T%ab(2,2)*T%ab(3,3)-T%ab(2,3)*T%ab(3,2))
-        inv_2%ab(2,1)=-1./detT*(T%ab(2,1)*T%ab(3,3)-T%ab(2,3)*T%ab(3,1))
-        inv_2%ab(3,1)=+1./detT*(T%ab(2,1)*T%ab(3,2)-T%ab(2,2)*T%ab(3,1))
-        inv_2%ab(1,2)=-1./detT*(T%ab(1,2)*T%ab(3,3)-T%ab(1,3)*T%ab(3,2))
-        inv_2%ab(2,2)=+1./detT*(T%ab(1,1)*T%ab(3,3)-T%ab(1,3)*T%ab(3,1))
-        inv_2%ab(3,2)=-1./detT*(T%ab(1,1)*T%ab(3,2)-T%ab(1,2)*T%ab(3,1))
-        inv_2%ab(1,3)=+1./detT*(T%ab(1,2)*T%ab(2,3)-T%ab(1,3)*T%ab(2,2))
-        inv_2%ab(2,3)=-1./detT*(T%ab(1,1)*T%ab(2,3)-T%ab(1,3)*T%ab(2,1))
-        inv_2%ab(3,3)=+1./detT*(T%ab(1,1)*T%ab(2,2)-T%ab(1,2)*T%ab(2,1))
+        inv_2%ab(1,1)=+idetT*(T%ab(2,2)*T%ab(3,3)-T%ab(2,3)*T%ab(3,2))
+        inv_2%ab(2,1)=-idetT*(T%ab(2,1)*T%ab(3,3)-T%ab(2,3)*T%ab(3,1))
+        inv_2%ab(3,1)=+idetT*(T%ab(2,1)*T%ab(3,2)-T%ab(2,2)*T%ab(3,1))
+        inv_2%ab(1,2)=-idetT*(T%ab(1,2)*T%ab(3,3)-T%ab(1,3)*T%ab(3,2))
+        inv_2%ab(2,2)=+idetT*(T%ab(1,1)*T%ab(3,3)-T%ab(1,3)*T%ab(3,1))
+        inv_2%ab(3,2)=-idetT*(T%ab(1,1)*T%ab(3,2)-T%ab(1,2)*T%ab(3,1))
+        inv_2%ab(1,3)=+idetT*(T%ab(1,2)*T%ab(2,3)-T%ab(1,3)*T%ab(2,2))
+        inv_2%ab(2,3)=-idetT*(T%ab(1,1)*T%ab(2,3)-T%ab(1,3)*T%ab(2,1))
+        inv_2%ab(3,3)=+idetT*(T%ab(1,1)*T%ab(2,2)-T%ab(1,2)*T%ab(2,1))
 
        end function inv_2
 
@@ -24,16 +24,16 @@ function inv_2s(T)
 
         type(Tensor2s), intent(in) :: T
         type(Tensor2s) :: inv_2s
-        real(kind=8) :: detT
+        real(kind=8) :: idetT
 
-        detT = det(T)
+        idetT = 1.d0/det(T)
 
-        inv_2s%a6(1)=+1./detT*(T%a6(2)*T%a6(3)  -T%a6(5)*T%a6(5))
-        inv_2s%a6(4)=-1./detT*(T%a6(4)*T%a6(3)  -T%a6(5)*T%a6(6))
-        inv_2s%a6(6)=+1./detT*(T%a6(4)*T%a6(5)  -T%a6(2)*T%a6(6))
-        inv_2s%a6(2)=+1./detT*(T%a6(1)*T%a6(3)  -T%a6(6)*T%a6(6))
-        inv_2s%a6(5)=-1./detT*(T%a6(1)*T%a6(5)  -T%a6(4)*T%a6(6))
-        inv_2s%a6(3)=+1./detT*(T%a6(1)*T%a6(2)  -T%a6(4)*T%a6(4))
+        inv_2s%a6(1)=+idetT*(T%a6(2)*T%a6(3)  -T%a6(5)*T%a6(5))
+        inv_2s%a6(4)=-idetT*(T%a6(4)*T%a6(3)  -T%a6(5)*T%a6(6))
+        inv_2s%a6(6)=+idetT*(T%a6(4)*T%a6(5)  -T%a6(2)*T%a6(6))
+        inv_2s%a6(2)=+idetT*(T%a6(1)*T%a6(3)  -T%a6(6)*T%a6(6))
+        inv_2s%a6(5)=-idetT*(T%a6(1)*T%a6(5)  -T%a6(4)*T%a6(6))
+        inv_2s%a6(3)=+idetT*(T%a6(1)*T%a6(2)  -T%a6(4)*T%a6(4))
 
        end function inv_2s
 
@@ -42,17 +42,19 @@ function inv2d(T,detT)
 
         type(Tensor2), intent(in) :: T
         type(Tensor2) :: inv2d
-        real(kind=8) :: detT
-
-        inv2d%ab(1,1)=+1./detT*(T%ab(2,2)*T%ab(3,3)-T%ab(2,3)*T%ab(3,2))
-        inv2d%ab(2,1)=-1./detT*(T%ab(2,1)*T%ab(3,3)-T%ab(2,3)*T%ab(3,1))
-        inv2d%ab(3,1)=+1./detT*(T%ab(2,1)*T%ab(3,2)-T%ab(2,2)*T%ab(3,1))
-        inv2d%ab(1,2)=-1./detT*(T%ab(1,2)*T%ab(3,3)-T%ab(1,3)*T%ab(3,2))
-        inv2d%ab(2,2)=+1./detT*(T%ab(1,1)*T%ab(3,3)-T%ab(1,3)*T%ab(3,1))
-        inv2d%ab(3,2)=-1./detT*(T%ab(1,1)*T%ab(3,2)-T%ab(1,2)*T%ab(3,1))
-        inv2d%ab(1,3)=+1./detT*(T%ab(1,2)*T%ab(2,3)-T%ab(1,3)*T%ab(2,2))
-        inv2d%ab(2,3)=-1./detT*(T%ab(1,1)*T%ab(2,3)-T%ab(1,3)*T%ab(2,1))
-        inv2d%ab(3,3)=+1./detT*(T%ab(1,1)*T%ab(2,2)-T%ab(1,2)*T%ab(2,1))
+        real(kind=8) :: detT, idetT
+
+        idetT = 1.d0/detT
+
+        inv2d%ab(1,1)=+idetT*(T%ab(2,2)*T%ab(3,3)-T%ab(2,3)*T%ab(3,2))
+        inv2d%ab(2,1)=-idetT*(T%ab(2,1)*T%ab(3,3)-T%ab(2,3)*T%ab(3,1))
+        inv2d%ab(3,1)=+idetT*(T%ab(2,1)*T%ab(3,2)-T%ab(2,2)*T%ab(3,1))
+        inv2d%ab(1,2)=-idetT*(T%ab(1,2)*T%ab(3,3)-T%ab(1,3)*T%ab(3,2))
+        inv2d%ab(2,2)=+idetT*(T%ab(1,1)*T%ab(3,3)-T%ab(1,3)*T%ab(3,1))
+        inv2d%ab(3,2)=-idetT*(T%ab(1,1)*T%ab(3,2)-T%ab(1,2)*T%ab(3,1))
+        inv2d%ab(1,3)=+idetT*(T%ab(1,2)*T%ab(2,3)-T%ab(1,3)*T%ab(2,2))
+        inv2d%ab(2,3)=-idetT*(T%ab(1,1)*T%ab(2,3)-T%ab(1,3)*T%ab(2,1))
+        inv2d%ab(3,3)=+idetT*(T%ab(1,1)*T%ab(2,2)-T%ab(1,2)*T%ab(2,1))
 
        end function inv2d
 
@@ -61,13 +63,15 @@ function inv2sd(T,detT)
 
         type(Tensor2s), intent(in) :: T
         type(Tensor2s) :: inv2sd
-        real(kind=8):: detT
-
-        inv2sd%a6(1)=+1./detT*(T%a6(2)*T%a6(3)  -T%a6(5)*T%a6(5))
-        inv2sd%a6(4)=-1./detT*(T%a6(4)*T%a6(3)  -T%a6(5)*T%a6(6))
-        inv2sd%a6(6)=+1./detT*(T%a6(4)*T%a6(5)  -T%a6(2)*T%a6(6))
-        inv2sd%a6(2)=+1./detT*(T%a6(1)*T%a6(3)  -T%a6(6)*T%a6(6))
-        inv2sd%a6(5)=-1./detT*(T%a6(1)*T%a6(5)  -T%a6(4)*T%a6(6))
-        inv2sd%a6(3)=+1./detT*(T%a6(1)*T%a6(2)  -T%a6(4)*T%a6(4))
+        real(kind=8):: detT, idetT
+
+        idetT = 1.d0/detT
+
+        inv2sd%a6(1)=+idetT*(T%a6(2)*T%a6(3)  -T%a6(5)*T%a6(5))
+        inv2sd%a6(4)=-idetT*(T%a6(4)*T%a6(3)  -T%a6(5)*T%a6(6))
+        inv2sd%a6(6)=+idetT*(T%a6(4)*T%a6(5)  -T%a6(2)*T%a6(6))
+        inv2sd%a6(2)=+idetT*(T%a6(1)*T%a6(3)  -T%a6(6)*T%a6(6))
+        inv2sd%a6(5)=-idetT*(T%a6(1)*T%a6(5)  -T%a6(4)*T%a6(6))
+        inv2sd%a6(3)=+idetT*(T%a6(1)*T%a6(2)  -T%a6(4)*T%a6(4))
 
        end function inv2sd