error DetectExtrudeStructure

Extension of Elmer in computational glaciology
Post Reply
nate maier
Posts: 7
Joined: 06 Jan 2017, 01:33
Antispam: Yes

error DetectExtrudeStructure

Post by nate maier »

Hey all -

I am doing some simple flow-line modelling in a relatively small model domain. I created an extruded mesh using gmsh. Any time I make a mesh with a resolution of less than 10 meters for either my z or x resolution, I get thrown the error DetectExtrudedStructure when I run my . sif file.

Code: Select all

MAIN: Reading Model: cavity_base.sif
LoadInputFile: Scanning input file: cavity_base.sif
LoadInputFile: Loading input file: cavity_base.sif
Model Input:  Unlisted keyword: [depth] in section: [initial condition 1]
Model Input:  Unlisted keyword: [mesh velocity first zero] in section: [solver 1]
Model Input:  Unlisted keyword: [flow solver name] in section: [solver 4]
Model Input:  Unlisted keyword: [strainrate variable name] in section: [solver 4]
Model Input:  Unlisted keyword: [rate factor 1] in section: [material 1]
Model Input:  Unlisted keyword: [rate factor 2] in section: [material 1]
Model Input:  Unlisted keyword: [activation energy 1] in section: [material 1]
Model Input:  Unlisted keyword: [activation energy 2] in section: [material 1]
Model Input:  Unlisted keyword: [limit temperature] in section: [material 1]
Model Input:  Unlisted keyword: [constant temperature] in section: [material 1]
Loading user function library: [StructuredMeshMapper]...[StructuredMeshMapper_Init0]
Loading user function library: [StructuredProjectToPlane]...[StructuredProjectToPlane_Init0]
Loading user function library: [ElmerIceSolvers]...[ComputeStrainRate_Init0]
LoadMesh: Base mesh name: ./ed_domain_2m
MAIN: -------------------------------------
AddVtuOutputSolverHack: Adding ResultOutputSolver to write VTU output in file: s40_a10_b120
Loading user function library: [StructuredMeshMapper]...[StructuredMeshMapper_Init]
Loading user function library: [StructuredMeshMapper]...[StructuredMeshMapper_bulk]
Loading user function library: [StructuredMeshMapper]...[StructuredMeshMapper]
Loading user function library: [StructuredProjectToPlane]...[StructuredProjectToPlane_Init]
Loading user function library: [StructuredProjectToPlane]...[StructuredProjectToPlane_bulk]
Loading user function library: [StructuredProjectToPlane]...[StructuredProjectToPlane]
Loading user function library: [FlowSolve]...[FlowSolver_Init]
Loading user function library: [FlowSolve]...[FlowSolver_bulk]
Loading user function library: [FlowSolve]...[FlowSolver]
Loading user function library: [ElmerIceSolvers]...[ComputeStrainRate_Init]
Loading user function library: [ElmerIceSolvers]...[ComputeStrainRate_bulk]
Loading user function library: [ElmerIceSolvers]...[ComputeStrainRate]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_Init]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver_bulk]
Loading user function library: [ResultOutputSolve]...[ResultOutputSolver]
StructuredMeshMapper: ---------------------------------------
StructuredMeshMapper: Performing mapping on a structured mesh
StructuredMeshMapper: ---------------------------------------
DetectExtrudedStructure:  Top and bottom pointer init time:   0.69200000000000017
DetectExtrudedStructure:  Top and bottom pointer init rounds:           90
DetectExtrudedStructure:  Number of nodes at the top:        10001
DetectExtrudedStructure:  Number of nodes at the bottom:        10001
StructuredMeshMapper:  Active coordinate mapping time:    3.2000000000000028E-002
MAIN: 
MAIN: -------------------------------------
MAIN:  Steady state iteration:            1
MAIN: -------------------------------------
MAIN: 
StructuredProjectToPlane: ------------------------------------------
StructuredProjectToPlane: Performing projection on a structured mesh
StructuredProjectToPlane: ------------------------------------------
DetectExtrudedStructure: Try to increase value for > Dot Product Tolerance <
ERROR:: DetectExtrudedStructure: Zero rounds implies unsuccesfull operation
I have modified my dot product tolerance - but can't seem to find any value that will make the DetectExtrudedStructure happy. I have attached my .sif file and my .geo file originally used to produce my mesh. Any help would be much appreciated.

Cheers,

Nate Maier

Code: Select all

!1/23/17 - cavity_base- 2-d flowline model to test how short length scale contrasts on in basal slipperiness will modify the basa!strain field. 


!echo on
Header
  !CHECK KEYWORDS Warn
  Mesh DB "." "ed_domain_2m" ! 10000 x 750 m model domain with 20 m resolution in x and y.
  Include Path ""
  Results Directory ""
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Simulation
  Max Output Level = 4
  Coordinate System = "Cartesian 2D"
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = "Steady"
  Steady State Max Iterations = 1
  Output Intervals = 1
  Output File = "s40_a10_b120.result" ! Need to edit the slip ratio in active region
  Post File = "s40_a10_b120.vtu" ! use .ep suffix for leagcy format
  Initialize Dirichlet Conditions = Logical False
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Constants
  Gas Constant = Real 8.314 !Joule/mol
  Stefan Boltzmann = 5.67e-08 ! I think for heat equation?
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!SLIP FUNCTION - this is our basal boundary slip function. Coded in fortran 90. 
$ function slide(X) {\
  if (X(0) > 5100)\
    {_slide = 0.00018;}\
  else if (X(0) < 4900)\
    {_slide = 0.00018;}\
  else\
   { _slide = 0.005;}\
}

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!Glens Flow Law - different way to apply glens flow law used below
$ function glen(Th) {\
   EF = 1.0;\
   AF = getArrheniusFactor(Th);\
   _glen = (2.0*EF*AF)^(-1.0/3.0);\
}
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!Arrhenius Factor Needed by Glens Flow Law - different way to apply glens flow law used below
$ function getArrheniusFactor(Th){ \
    if (Th<-10) {_getArrheniusFactor=3.985E-13 * exp( -60.0E03/(8.314 * (273.15 + Th)));}\
    else {\
       if (Th>0) _getArrheniusFactor=1.916E03 * exp( -139.0E03/(8.314 *  (273.15)));\
            else _getArrheniusFactor=1.916E03 * exp( -139.0E03/(8.314 *  (273.15 + Th)));}\
}
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Body 1
  Name = "Glacier"
  Body Force = 1
  Equation = 1
  Material = 1
  Initial Condition = 1
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Equation 1
  Name = "Equation1"
  Convection = "computed"
  Flow Solution Name = String "Flow Solution"
  Active Solvers(4) = 1 2 3 4
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Initial Condition 1
!  Velocity 1 = 0.0
!  Velocity 2 = 0.0
!  Pressure = 0.0
  Depth = Real 0.0
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Is used to map surface and basal boundaries into mesh
Solver 1
  Exec Solver = "before Simulation"
  Equation = "MapCoordinate"
  Procedure = "StructuredMeshMapper" "StructuredMeshMapper"
  Active Coordinate = Integer 2 ! the mesh-update is y-direction
! For time being this is currently externally allocated
  Mesh Velocity Variable = String "Mesh Velocity 2"
! The 1st value is special as the mesh velocity could be unrelistically high
  Mesh Velocity First Zero = Logical True
! The accuracy applied to vector-projections
  Dot Product Tolerance = Real 0.01
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Solves for height above the bed and depth below the free surface
Solver 2
  Equation = "HeightDepth"
  Procedure = "StructuredProjectToPlane" "StructuredProjectToPlane"
  Active Coordinate = Integer 2
  Operator 1 = depth
  Operator 2 = height
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! the central part of the problem: the Stokes solver
Solver 3
!  Exec Solver = "Never" # uncommenting would switch this off 
  Equation = "Navier-Stokes"
  Optimize Bandwidth = Logical True
  ! direct solver
  Linear System Solver = Direct
  Linear System Direct Method = "UMFPACK"     
  ! alternative to above - Krylov subspace iterative solution
!  Linear System Solver = "Iterative"
!  Linear System Iterative Method =  "GCR"     !or "BICGStab"
  Linear System Max Iterations = 5000
  Linear System Convergence Tolerance = 1.0E-06
  Linear System Abort Not Converged = False
  Linear System Preconditioning = "ILU1"
  Linear System Residual Output = 1
  Flow Model = Stokes

  Steady State Convergence Tolerance = 1.0E-05
!  Stabilization Method can be [Stabilized,P2/P1,Bubbles] 
  Stabilization Method = Bubbles

  Nonlinear System Convergence Tolerance = 5.0E-04 
  Nonlinear System Convergence Measure = Solution
  Nonlinear System Max Iterations = 100
  Nonlinear System Newton After Iterations = 10
  Nonlinear System Newton After Tolerance =  1.0E-01
  Nonlinear System Relaxation Factor = 0.75

  Exported Variable 1 = -dofs 3 "Mesh Velocity" 
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!NM - Strain rate solver - exports the strain-rate tensor in results file
Solver 4
  Equation = "Eij"          
  Variable = -nooutput "Eij"     
  Variable DOFs = 1

  Exported Variable 1 = "StrainRate"
  Exported Variable 1 DOFs = 5

  Procedure = "ElmerIceSolvers" "ComputeStrainRate"
  Flow Solver Name = String "Flow Solution"
  StrainRate Variable Name = String "StrainRate"

  Linear System Solver = Direct         
  Linear System Direct Method = umfpack
End 

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! we use m-yr-MPa system 1 yr = 31556926.0 sec
! need to spend some time looking at this
Material 1
  Name = "ice"
  Density = Real $910.0*1.0E-06*(31556926.0)^(-2.0)
  Viscosity Model = String "Glen"
  ! Glen's flow law (using Glen)
!----------------
! viscosity stuff
!----------------
  Viscosity Model = String "Glen"
! Viscosity has to be set to a dummy value
! to avoid warning output from Elmer
  Viscosity = Real 1.0 
  Glen Exponent = Real 3.0
  Critical Shear Rate = Real 1.0e-10
! Rate factors (Paterson value in MPa^-3a^-1)
  Rate Factor 1 = Real 1.258e13  
  Rate Factor 2 = Real 6.046e28
! these are in SI units - no problem, as long as
! the gas constant also is 
  Activation Energy 1 = Real 60e3
  Activation Energy 2 = Real 139e3  
  Glen Enhancement Factor = Real 1.0
! the variable taken to evaluate the Arrhenius law
! in general this should be the temperature relative
! to pressure melting point. The suggestion below plugs
! in the correct value obtained with TemperateIceSolver
  Temperature Field Variable = String "Temp Homologous"
! the temperature to switch between the 
! two regimes in the flow law
  Limit Temperature = Real -10.0
!
  Constant Temperature = Real -3.0
! Setting temperature to vary as a function of height using data
!!NOTE!! The lines below have to be justified exactly like this or they will not work
  !Constant Temperature = Variable height
  !Real cubic
     !include "GL15SiteCA_Temp_form.dat"
  !End
End
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!Defines gravity as body force
Body Force 1
  Name = "BodyForce1"
  Heat Source = 1
  Flow BodyForce 1 = Real 0.0                          
  Flow BodyForce 2 = Real $-9.81 * (31556926.0)^(2.0)  !MPa - a - m
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Boundary Condition 1
  Name = "bedrock"
  Target Boundaries = 1
! Using MATC to set bedrock boundry to bottom of mesh
  Bottom Surface = Variable Coordinate 1  
     Real MATC "200-tx*tan(0.0174533)+2.5*sin((2*3.14/20)*tx)"
  Normal-Tangential Velocity = Logical True
  Velocity 1 = Real 0.0e0   
  Slip Coefficient 2 = Real 0.0011
  !Slip Coefficient 2 = Variable Coordinate 1
     !include "perc_50.dat"
  !End
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Boundary Condition 2
  Name = "surface"
  Target Boundaries = 2
!Using MACT function to set surface slope based on x coordinate
  Top Surface = Variable Coordinate 1  
     Real MATC "900-tx*tan(0.0174533)"
End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!par!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Upgradient Boundary Condition - Perscribed Velocity - this boundary condition has trouble with convergence without loosening to tolerances for the stokes solver
Boundary Condition 3
  Name = "sides"
  Target Boundaries(2) = 3 4 ! combine left and right boundary
!Velocity 1 = Real 300
  Velocity 1 = Variable Depth
  Real cubic
      include "bound_120.dat"
      


!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

! Downgradient Boundary Condition - External Pressure
!Boundary Condition 4
! Name = "downgradient"
  !Target Boundaries = 4 
  !Velocity 1 = Real 300
  !Velocity 1 = Variable Depth
  !Real cubic
     !include "bound_down_120.dat"
  !End
!End

Code: Select all

Point(1) = {10000, 0, 0, 2};
Point(2) = {0, 0, 0, 2};
Line(1) = {2, 1};
Extrude {0, 900, 0} {
  Line{1};Layers{450};Recombine;
}
Physical Surface(6) = {5};
Physical Line(7) = {1};
Physical Line(8) = {2};
Physical Line(9) = {3};
Physical Line(10) = {4};
tzwinger
Site Admin
Posts: 99
Joined: 24 Aug 2009, 12:20
Antispam: Yes

Re: error DetectExtrudeStructure

Post by tzwinger »

Hi,
given the fact that it works with coarser resolutions, I'd try to provide the

Code: Select all

Dot Product Tolerance
keyword also in the "StructuredToPlane" solver section.

Best Regards,

Thomas
nate maier
Posts: 7
Joined: 06 Jan 2017, 01:33
Antispam: Yes

Re: error DetectExtrudeStructure

Post by nate maier »

Thomas -

I tried putting the keyword into the StructuretoPlane solver section. If I increase the dot product tolerance to values > 10 the .sif runs, but only goes through one iteration and outputs garbage. That would be expected, poking around elsewhere in the elmer forums, the dot product tolerance should be <<1. Any other recommendations?

Cheers,

Nate
raback
Site Admin
Posts: 4812
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: error DetectExtrudeStructure

Post by raback »

Hi

The default is 1.0e-4.

The code tests that an edge is oriented in the direction of the extrusion. Then the projection from the dot product is one. It is tolerance to this 1 that is checked for. It does not make sense to go to ~1 or so. That basically would mean that all nodes within element are aligned.

Why the tolerance changes may be related to number of digits. Assume that you have rounding error to different digit. With small mesh size parameter h vs. the relative rounding error gets bigger and the test is failed. So you could try to increase number of digits also. In ElmerGrid you do this with "-decimals" flag. It could be also something already in the gmsh mesh.

-Peter
tzwinger
Site Admin
Posts: 99
Joined: 24 Aug 2009, 12:20
Antispam: Yes

Re: error DetectExtrudeStructure

Post by tzwinger »

Hi,

to test Peter's suspicion, whether the problem is from gmsh, you could easily create this simple geometry using native format of ElmerGrid. To find a blueprint, just check for any square.grd in the existing Elmer examples.

Regards,

Thomas
nate maier
Posts: 7
Joined: 06 Jan 2017, 01:33
Antispam: Yes

Re: error DetectExtrudeStructure

Post by nate maier »

Let me do that and get back to you both.

Thanks,

Nate
nate maier
Posts: 7
Joined: 06 Jan 2017, 01:33
Antispam: Yes

Re: error DetectExtrudeStructure

Post by nate maier »

Hey all -

I tried adding the -decimal flag and adding more decimal places when converting my gmsh generated mesh, .msh to the .mesh for running with the .sif file. No luck there.

I also tried setting the dot product tolerance to the default value of 1e-4 for both the StructuredMeshMapper and StructuredProjectoPlane subsections of my .sif file. No luck there either.

I then generated a mesh in the native .grd format for ElmerGrid. Right now it seems to be working quite well. I no longer get the error DetectExtrudeStructure when running my .sif file for a high resolution mesh. I will probably inspect my .mesh files generated from the gmsh. I am curious as to the cause, in the long runI think gmsh is a bit more capable then ElmerGrid for mesh generation so it would be worth while figuring out the issue.

Thanks so much for the help Thomas and Peter,

Nate
raback
Site Admin
Posts: 4812
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: error DetectExtrudeStructure

Post by raback »

Hi

Probably the initial triangular mesh is not perfectly extruded to the last decimal. You could have tried sloppier tolerance, for example 1.0e-2, but not too sloppy like 1.

-Peter
Post Reply