"unitlength" parameter, which determines the physical size of the mesh, now available as a global mesh variable "mesh_unitlength" and written to the .mesh file during init.

Hence, it is available to the marc_deformedGeometry.py script via "mesh_init_postprocessing()" and "mesh_get_unitlength()", so no need for setting the scaling of the displacement vectors explicitly through an option; now displacements and nodal positions are always consistent

code/damask.core.pyf

@@ -203,6 +203,10 @@ python module core ! in
             integer :: mesh_get_Ncellnodes
         end function mesh_get_Ncellnodes
         
+        function mesh_get_unitlength() ! in :mesh:mesh.f90
+            real*8 :: mesh_get_unitlength
+        end function mesh_get_unitlength
+        
         function mesh_get_nodeAtIP(elemtypeFE,ip) ! in :mesh:mesh.f90
             character(len=*), intent(in) :: elemtypeFE
             integer,          intent(in) :: ip

code/mesh.f90

@@ -57,6 +57,9 @@ module mesh
  integer(pInt), dimension(:,:,:,:), allocatable, public, protected :: &
    mesh_ipNeighborhood                                                                              !< 6 or less neighboring IPs as [element_num, IP_index, neighbor_index that points to me]
 
+ real(pReal), public, protected :: &
+   mesh_unitlength                                                                                  !< physical length of one unit in mesh
+
  real(pReal), dimension(:,:), allocatable, public :: &
    mesh_node, &                                                                                     !< node x,y,z coordinates (after deformation! ONLY FOR MARC!!!)
    mesh_cellnode                                                                                    !< cell node x,y,z coordinates (after deformation! ONLY FOR MARC!!!)
@@ -421,6 +424,7 @@ module mesh
    mesh_cellCenterCoordinates, &
    mesh_init_postprocessing, &
    mesh_get_Ncellnodes, &
+   mesh_get_unitlength, &
    mesh_get_nodeAtIP
 #ifdef Spectral
  public :: &
@@ -512,7 +516,8 @@ subroutine mesh_init(ip,el)
    IO_open_InputFile
 #endif
  use numerics, only: &
-   usePingPong
+   usePingPong, &
+   numerics_unitlength
  use FEsolving, only: &
    FEsolving_execElem, &
    FEsolving_execIP, &
@@ -549,7 +554,10 @@ subroutine mesh_init(ip,el)
  if (allocated(FE_ipNeighbor))                deallocate(FE_ipNeighbor)
  if (allocated(FE_cellnodeParentnodeWeights)) deallocate(FE_cellnodeParentnodeWeights)
  if (allocated(FE_subNodeOnIPFace))           deallocate(FE_subNodeOnIPFace)
+
  call mesh_build_FEdata                                                                             ! get properties of the different types of elements
+ mesh_unitlength = numerics_unitlength                                                              ! set physical extent of a length unit in mesh
+
 #ifdef Spectral
  call IO_open_file(fileUnit,geometryFile)                                                           ! parse info from geometry file...
  res = mesh_spectral_getGrid(fileUnit)
@@ -1265,8 +1273,6 @@ end subroutine mesh_spectral_count_cpSizes
 !! Allocates global arrays 'mesh_node0' and 'mesh_node'
 !--------------------------------------------------------------------------------------------------
 subroutine mesh_spectral_build_nodes()
- use numerics, &
-   only: numerics_unitlength
 
  implicit none
  integer(pInt) :: n
@@ -1275,13 +1281,13 @@ subroutine mesh_spectral_build_nodes()
  allocate ( mesh_node  (3,mesh_Nnodes) ); mesh_node  = 0.0_pReal
  
  forall (n = 0_pInt:mesh_Nnodes-1_pInt)
-   mesh_node0(1,n+1_pInt) = numerics_unitlength * &
+   mesh_node0(1,n+1_pInt) = mesh_unitlength * &
            geomdim(1) * real(mod(n,(res(1)+1_pInt) ),pReal) &
                                                    / real(res(1),pReal)
-   mesh_node0(2,n+1_pInt) = numerics_unitlength * &
+   mesh_node0(2,n+1_pInt) = mesh_unitlength * &
            geomdim(2) * real(mod(n/(res(1)+1_pInt),(res(2)+1_pInt)),pReal) &
                                                    / real(res(2),pReal)
-   mesh_node0(3,n+1_pInt) = numerics_unitlength * &
+   mesh_node0(3,n+1_pInt) = mesh_unitlength * &
            geomdim(3) * real(mod(n/(res(1)+1_pInt)/(res(2)+1_pInt),(res(3)+1_pInt)),pReal) &
                                                    / real(res(3),pReal)
  end forall 
@@ -2762,7 +2768,6 @@ subroutine mesh_marc_build_nodes(myUnit)
                  IO_stringPos, &
                  IO_fixedIntValue, &
                  IO_fixedNoEFloatValue
-use numerics, only: numerics_unitlength
 
  implicit none
  integer(pInt), intent(in) :: myUnit
@@ -2788,7 +2793,7 @@ use numerics, only: numerics_unitlength
        read (myUnit,610,END=670) line
        m = mesh_FEasCP('node',IO_fixedIntValue(line,node_ends,1_pInt))
        do j = 1_pInt,3_pInt
-         mesh_node0(j,m) = numerics_unitlength * IO_fixedNoEFloatValue(line,node_ends,j+1_pInt)
+         mesh_node0(j,m) = mesh_unitlength * IO_fixedNoEFloatValue(line,node_ends,j+1_pInt)
        enddo
      enddo
      exit
@@ -3419,7 +3424,6 @@ subroutine mesh_abaqus_build_nodes(myUnit)
                  IO_error, &
                  IO_countDataLines, &
                  IO_intValue
-use numerics, only: numerics_unitlength
 
  implicit none
  integer(pInt), intent(in) :: myUnit
@@ -3460,7 +3464,7 @@ use numerics, only: numerics_unitlength
        myPos = IO_stringPos(line,maxNchunks)
        m = mesh_FEasCP('node',IO_intValue(line,myPos,1_pInt))
        do j=1_pInt, 3_pInt
-         mesh_node0(j,m) = numerics_unitlength * IO_floatValue(line,myPos,j+1_pInt)
+         mesh_node0(j,m) = mesh_unitlength * IO_floatValue(line,myPos,j+1_pInt)
        enddo  
      enddo
    endif
@@ -5212,6 +5216,7 @@ subroutine mesh_write_meshfile
  integer(pInt) :: e,i,t,g,c,n
 
  call IO_write_jobFile(fileUnit,'mesh') 
+ write(fileUnit,'(A16,E10)') 'unitlength', mesh_unitlength
  write(fileUnit,'(A16,I10)') 'maxNcellnodes', mesh_maxNcellnodes
  write(fileUnit,'(A16,I10)') 'maxNips', mesh_maxNips
  write(fileUnit,'(A16,I10)') 'maxNnodes', mesh_maxNnodes
@@ -5251,6 +5256,7 @@ integer function mesh_read_meshfile(filepath)
  integer(pInt) :: e,i,t,g,n
 
  open(fileUnit,status='old',err=100,iostat=mesh_read_meshfile,action='read',file=filepath)
+ read(fileUnit,'(TR16,E10)',err=100,iostat=mesh_read_meshfile) mesh_unitlength
  read(fileUnit,'(TR16,I10)',err=100,iostat=mesh_read_meshfile) mesh_maxNcellnodes
  read(fileUnit,'(TR16,I10)',err=100,iostat=mesh_read_meshfile) mesh_maxNips
  read(fileUnit,'(TR16,I10)',err=100,iostat=mesh_read_meshfile) mesh_maxNnodes
@@ -5310,6 +5316,18 @@ integer(pInt) function mesh_get_Ncellnodes()
 end function mesh_get_Ncellnodes
 
 
+!--------------------------------------------------------------------------------------------------
+!> @brief returns global variable mesh_unitlength
+!--------------------------------------------------------------------------------------------------
+real(pReal) function mesh_get_unitlength()
+
+ implicit none
+ 
+ mesh_get_unitlength = mesh_unitlength
+
+end function mesh_get_unitlength
+
+
 !--------------------------------------------------------------------------------------------------
 !> @brief returns node that is located at an ip
 !> @details return zero if requested ip does not exist or not available (more ips than nodes)

processing/post/marc_deformedGeometry.py

@@ -42,13 +42,10 @@ parser.add_option('-r','--range', dest='range', type='int', nargs=3, \
                   help='range of positions (or increments) to output (start, end, step) [all]')
 parser.add_option('--increments', action='store_true', dest='getIncrements', \
                   help='switch to increment range [%default]')
-parser.add_option('--scale', dest='scale', type='float', \
-                  help='scaling factor for the nodal displacements [%default]')
 
 
 parser.set_defaults(dir = 'vtk')
 parser.set_defaults(getIncrements= False)
-parser.set_defaults(scale = 1)
 
 (options, files) = parser.parse_args()
 
@@ -84,6 +81,7 @@ if damask.core.mesh.mesh_init_postprocessing(filename+'.mesh') > 0:
   print('error: init not successful')
   sys.exit(-1)
 Ncellnodes = damask.core.mesh.mesh_get_Ncellnodes()
+unitlength = damask.core.mesh.mesh_get_unitlength()
 
 
 # ---------------------------   create output dir   --------------------------------
@@ -131,7 +129,7 @@ for incCount,position in enumerate(locations):     # walk through locations
   node_displacement = [[0,0,0] for i in range(Nnodes)]
   for n in range(Nnodes):
     if p.node_displacements():
-      node_displacement[n] = map(lambda x:x*options.scale,list(p.node_displacement(n)))
+      node_displacement[n] = map(lambda x:x*unitlength,list(p.node_displacement(n)))
   c = damask.core.mesh.mesh_build_cellnodes(numpy.array(node_displacement).T,Ncellnodes)
   cellnode_displacement = [[c[i][n] for i in range(3)] for n in range(Ncellnodes)]