Skip to content

Commit

Permalink
Merge branch 'devel' of https://github.com/ElmerCSC/elmerfem into devel
Browse files Browse the repository at this point in the history
  • Loading branch information
ettaka committed Nov 15, 2019
2 parents 60b5880 + 4b4b4de commit 1bd8687
Show file tree
Hide file tree
Showing 46 changed files with 2,510 additions and 717 deletions.
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
sif/cylinders.sif
Original file line number Diff line number Diff line change
@@ -0,0 +1,75 @@
# A test case for FreeCAD automatic scripting with thermo-structural (Elmer test is a modified fem/tests/ThermalBiMetal2)
# original date: November 2019
# Author: Eelis Takala
# email: eelis.takala@gmail.com
doc = App.newDocument('Cylinders')
import sys
import os
import time
import math
import FreeCAD
import Part


sys.path.insert(0, os.path.join(os.path.dirname(__file__), '..', '..'))
import FreeCADBatchFEMTools as FBFT

script_directory = os.path.dirname(__file__)

def cylinder(r1, r2, h, name='cyl', mesh_size=10.):

cyl = doc.addObject("PartDesign::Body","Cylinder")
sketch = cyl.newObject('Sketcher::SketchObject',name+' sketch')
sketch.MapMode = 'FlatFace'
sketch.addGeometry(Part.Circle(App.Vector(0.000000,0.000000,0),App.Vector(0,0,1),r2),False)
sketch.addGeometry(Part.Circle(App.Vector(0.000000,0.000000,0),App.Vector(0,0,1),r1),False)
doc.recompute()

pad = cyl.newObject("PartDesign::Pad",name + " pad")
pad.Profile = sketch
pad.Length = h
pad.Length2 = h
pad.Type = 0
pad.UpToFace = None
pad.Reversed = 0
pad.Midplane = 1
pad.Offset = 0.000000
doc.recompute()

# create entities dict
face_picks = [('b1', 0), ('b2', 1), ('b3', 2)]
faces = FBFT.pick_faces_from_geometry(cyl, face_picks)
solids = []
FBFT.add_entity_in_list(solids, name, cyl, {'mesh size': mesh_size})
entities_dict = FBFT.create_entities_dict(name, faces, solids, cyl)

return entities_dict

def create_geometry():

mesh_size_max = 100.
cyl_entities_list = []
cyl_entities_list.append(cylinder(20, 10, 40., name='cyl1', mesh_size=3.))
cyl_entities_list.append(cylinder(30, 20, 20., name='cyl2', mesh_size=4.))
entities_dict = FBFT.merge_entities_dicts(cyl_entities_list, 'All', default_mesh_size=mesh_size_max,
add_prefixes={'solids': False, 'faces': True})
solid_objects = FBFT.get_solids_from_entities_dict(entities_dict)

mesh_object, compound_filter = FBFT.create_mesh_object_and_compound_filter(solid_objects, mesh_size_max, doc)

FBFT.find_boundaries_with_entities_dict(mesh_object, compound_filter, entities_dict, doc)
body_mesh_groups = FBFT.find_bodies_with_entities_dict(mesh_object, compound_filter, entities_dict, doc, point_search=False)
FBFT.define_mesh_sizes(mesh_object, compound_filter, entities_dict, doc, point_search=False, ignore_list=['air'])

FBFT.fit_view()
FBFT.create_mesh(mesh_object)
FBFT.export_unv(os.path.join(script_directory, 'cylinders.unv'), mesh_object)

try:
create_geometry()

except Exception:
import traceback
print(str(traceback.format_exc()))

exit()
Original file line number Diff line number Diff line change
@@ -0,0 +1,10 @@
#!/bin/bash

# A test case for FreeCAD automatic scripting with thermo-structural (Elmer test is a modified fem/tests/ThermalBiMetal2)
# original date: November 2019
# Author: Eelis Takala
# email: eelis.takala@gmail.com
FreeCAD -c ${PWD}/cylinders.py
ElmerGrid 8 2 cylinders.unv -autoclean
ElmerSolver
rm cylinders TEST.PASSED -r
Original file line number Diff line number Diff line change
@@ -0,0 +1,106 @@
! A test case for FreeCAD automatic scripting with thermo-structural (Elmer test is a modified fem/tests/ThermalBiMetal2)
! original date: November 2019
! Author: Eelis Takala
! email: eelis.takala@gmail.com
Check Keywords Warn
INCLUDE "cylinders/mesh.names"
Header
Mesh DB "." "cylinders"
End
Simulation
Max Output Level = 5
Coordinate System = Cartesian
Simulation Type = Steady State
Steady State Max Iterations = 1
! Post File = "case.ep"
Post File = cylinders.vtu
End
Constants
Permittivity Of Vacuum = 8.8542e-12
End
Body 1
Name = cyl1
Target Bodies(1) = $ cyl1
Equation = 1
Material = 1
Initial Condition = 1
End
Body 2
Name = cyl2
Target Bodies(1) = $ cyl2
Equation = 1
Material = 2
Body Force = 1
End
Body Force 1
Heat Source = 1.0e3
End
Equation 1
Active Solvers(2) = 1 2
Plane Stress = FALSE
End
Solver 1
Equation = HeatSolver
Variable = Temperature
Procedure = "HeatSolve" "HeatSolver"
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-4
Nonlinear System Max Iterations = 1
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-5
Nonlinear System Relaxation Factor = 1
Nonlinear System Convergence Measure = solution
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-8
Linear System Preconditioning = ILU1
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 10
Linear System Precondition Recompute = 1
End
Solver 2
Equation = "LinearDisp"
Procedure = "StressSolve" "StressSolver"
Variable = "Displacement"
Variable DOFs = Integer 3
Linear System Solver = Direct
Linear System Symmetric = Logical True
Linear System Scaling = Logical False
Linear System Iterative Method = BiCGStab
Linear System Direct Method = UMFPACK
Linear System Convergence Tolerance = 1.0e-8
Linear System Max Iterations = 200
Linear System Preconditioning = ILU2
Nonlinear System Convergence Tolerance = Real 1.0e-7
Nonlinear System Max Iterations = Integer 1
Nonlinear System Relaxation Factor = Real 1
Steady State Convergence Tolerance= 1.0e-6
Optimize Bandwidth = True
End
Material 1
Density = Real 1
Youngs Modulus = 1e9
Poisson Ratio = Real 0.3
Reference Temperature = 300.0
Heat Expansion Coefficient = 1.0e-4
Heat Conductivity = 1.0
End
Material 2
Density = Real 1
Youngs Modulus = 1e9
Poisson Ratio = Real 0.3
Reference Temperature = 300.0
Heat Expansion Coefficient = 2.0e-4
Heat Conductivity = 10.0
End
Boundary Condition 1
Name = cyl1_b1
Target Boundaries(1) = $ cyl1_b1
Displacement 1 = 0.0
Displacement 2 = 0.0
Temperature = 300.0
End
Solver 1 :: Reference Norm = 1.016E+04
Solver 1 :: Reference Norm Tolerance = 1.0E-3
10 changes: 5 additions & 5 deletions elmerice/permafrostmaterialdb.dat
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ Sand 1
0.35 eta0
0.1 etak
3.5 ks0th
0.006 ew
0.006 e1
0.0 bs
2600.0 rhos0
840.0 cs0
Expand Down Expand Up @@ -37,7 +37,7 @@ Granite 2
0.01 eta0
0.01 etak
3.0851 ks0th
0.006 ew
0.006 e1
0.6411 bs
2780.0 rhos0
700.0 cs0
Expand Down Expand Up @@ -70,7 +70,7 @@ Till 3
0.25 eta0
0.25 etak
2.0 ks0th
0.002 ew
0.002 e1
0.0 bs
2000.0 rhos0
1250.0 cs0
Expand Down Expand Up @@ -103,7 +103,7 @@ Gravel 4
0.25 eta0
0.25 etak
5.5 ks0th
0.006 ew
0.006 e1
0.0 bs
2650.0 rhos0
800.0 cs0
Expand Down Expand Up @@ -136,7 +136,7 @@ Silt 5
0.40 eta0
0.40 etak
1.8 ks0th
0.0025 ew
0.0025 e1
0.0 bs
2650.0 rhos0
900.0 cs0
Expand Down
43 changes: 25 additions & 18 deletions fem/src/BlockSolve.F90
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ MODULE BlockSolve
LOGICAL, PRIVATE :: isParallel=.FALSE.

TYPE(Variable_t), POINTER :: SolverVar => Null()
TYPE(Matrix_t), POINTER :: SolverMatrix => Null()
TYPE(Matrix_t), POINTER :: SolverMatrix => Null(), SaveMatrix

CONTAINS

Expand Down Expand Up @@ -582,11 +582,11 @@ SUBROUTINE BlockPickMatrix( Solver, NoVar )
IF( Amat % NumberOfRows > 0 ) THEN
SumAbsMat = SUM( ABS( Amat % Values ) )
IF( SumAbsMat < SQRT( TINY( SumAbsMat ) ) ) THEN
CALL Info('BlockSolver','Matrix is actually all zero, eliminating it!',Level=20)
CALL Info('BlockSolver','Matrix is actually all zero, eliminating it!',Level=12)
DEALLOCATE( Amat % Values )
IF( .NOT. ReuseMatrix ) THEN
DEALLOCATE( Amat % Rows, Amat % Cols )
IF( RowVar == ColVar ) DEALLOCATE( Amat % Diag )
IF( RowVar == ColVar ) DEALLOCATE( Amat % Diag, Amat % rhs )
END IF
Amat % NumberOfRows = 0
END IF
Expand Down Expand Up @@ -2346,7 +2346,7 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )
TYPE(Variable_t), POINTER :: Var, Var_save

LOGICAL :: GotOrder, BlockGS, Found, NS, ScaleSystem, DoSum, &
IsComplex, BlockScaling, DiagScaling, UsePrecMat
IsComplex, BlockScaling, DiagScaling, ThisScaling, UsePrecMat
CHARACTER(LEN=MAX_NAME_LEN) :: str
#ifndef USE_ISO_C_BINDINGS
INTEGER(KIND=AddrInt) :: AddrFunc
Expand Down Expand Up @@ -2393,8 +2393,9 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )

#define SOLSYS
#ifdef SOLSYS
IF (isParallel) &
IF (isParallel) THEN
ALLOCATE( x(TotMatrix % MaxSize), b(TotMatrix % MaxSize) )
END IF
#endif

! Initial guess
Expand All @@ -2416,8 +2417,10 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )
END IF

WRITE(Message,'(A,I0)') 'Solving block: ',i
CALL Info('BlockMatrixPrec',Message,Level=6)
CALL Info('BlockMatrixPrec',Message,Level=8)

CALL ListPushNameSpace('block '//TRIM(i2s(i))//TRIM(i2s(i))//':')

! Set pointers to the new linear system
!-------------------------------------------------------------------
Var => TotMatrix % SubVector(i) % Var
Expand Down Expand Up @@ -2481,8 +2484,6 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )
END DO
END IF

CALL ListPushNameSpace('block '//TRIM(i2s(i))//TRIM(i2s(i))//':')

! We do probably not want to compute the change within each iteration
CALL ListAddLogical( Asolver % Values,'Skip Advance Nonlinear iter',.TRUE.)
CALL ListAddLogical( Asolver % Values,'Skip Compute Nonlinear Change',.TRUE.)
Expand All @@ -2496,21 +2497,24 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )

IF( BlockScaling ) CALL BlockMatrixScaling(.TRUE.,i,i,b,UsePrecMat)

IF( DiagScaling .AND. UsePrecMat ) THEN
ThisScaling = ListGetLogical( Params,'Linear System Scaling',Found )
IF( .NOT. Found ) ThisScaling = DiagScaling

IF( ThisScaling .AND. UsePrecMat ) THEN
n = A % NumberOfRows
ALLOCATE( diagtmp(n), btmp(n) )

IF( TotMatrix % GotBlockStruct ) THEN
k = TotMatrix % InvBlockStruct(i)
IF( k <= 0 ) THEN
CALL Fatal('BlockMatrixPrec','Cannot define the originating block for scaling!')
CALL Fatal('BlockMatrixPrec','Cannot define the originating block '&
//TRIM(I2S(i))//' for scaling!')
END IF
l = SIZE( TotMatrix % BlockStruct )
ELSE
k = i
l = NoVar
END IF

l = Solver % Variable % DOFs
Diagtmp(1:n) = Solver % Matrix % DiagScaling(k::l)

! Scale x & b to the unscaled system of the tailored preconditioning matrix for given block.
Expand Down Expand Up @@ -2555,7 +2559,7 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )

IF( BlockScaling ) CALL BlockMatrixScaling(.FALSE.,i,i,b,UsePrecMat)

IF( DiagScaling .AND. UsePrecMat ) THEN
IF( ThisScaling .AND. UsePrecMat ) THEN
x(1:n) = x(1:n) / diagtmp(1:n)
DEALLOCATE( diagtmp, btmp )
END IF
Expand Down Expand Up @@ -2686,10 +2690,10 @@ SUBROUTINE BlockMatrixPrec( u,v,ipar )
Solver % Variable => Var_save

IF( BlockGS ) THEN
DEALLOCATE( vtmp, rtmp )
DEALLOCATE( vtmp, rtmp, xtmp )
END IF

CALL Info('BlockMatrixPrec','Finished block matrix preconditioning',Level=6)
CALL Info('BlockMatrixPrec','Finished block matrix preconditioning',Level=8)

END SUBROUTINE BlockMatrixPrec

Expand Down Expand Up @@ -2889,7 +2893,8 @@ SUBROUTINE BlockKrylovIter( Solver, MaxChange )
INTEGER :: NoVar, ndim, maxsize
LOGICAL :: Converged, Diverged
INTEGER :: Rounds, OutputInterval, PolynomialDegree
INTEGER, POINTER :: Offset(:),poffset(:),BlockStruct(:)
INTEGER, POINTER :: Offset(:),BlockStruct(:)
INTEGER, ALLOCATABLE, TARGET :: poffset(:)
INTEGER :: i,j,k,l,ia,ib,istat
LOGICAL :: LS, BlockAV,Found

Expand All @@ -2905,7 +2910,7 @@ SUBROUTINE BlockKrylovIter( Solver, MaxChange )
NoVar = TotMatrix % NoVar

CALL Info('BlockKrylovIter','Allocating temporal vectors for block system of size: '&
//TRIM(I2S(ndim)),Level=8)
//TRIM(I2S(ndim)),Level=10)

ALLOCATE(x(ndim), b(ndim),r(ndim),STAT=istat)
IF( istat /= 0 ) THEN
Expand Down Expand Up @@ -3136,7 +3141,7 @@ SUBROUTINE BlockSolveInt(A,x,b,Solver)
LOGICAL :: GotSlaveSolvers, SkipVar


TYPE(Matrix_t), POINTER :: Amat, SaveMatrix, SaveCM
TYPE(Matrix_t), POINTER :: Amat, SaveCM
TYPE(Mesh_t), POINTER :: Mesh
TYPE(ValueList_t), POINTER :: Params

Expand Down Expand Up @@ -3292,6 +3297,8 @@ SUBROUTINE BlockSolveInt(A,x,b,Solver)
ParEnv = Amat % ParMatrix % ParEnv
END IF
END IF

Solver % Variable => SolverVar
END DO
END DO
END IF
Expand Down
Loading

0 comments on commit 1bd8687

Please sign in to comment.