@@ -59,17 +59,17 @@ class CustomParser(
5959physical_parameters .add_argument ("--E" , dest = "E" , type = float , default = 2.0e4 , help = "Young's modulus" )
6060physical_parameters .add_argument ("--nu" , dest = "nu" , type = float , default = 0.3 , help = "Poisson's ratio" )
6161physical_parameters .add_argument (
62- "--disp" , type = float , default = - 0.2 , help = "Displacement in the y/z direction (2D/3D)"
62+ "--disp" , type = float , default = - 0.25 , help = "Displacement in the y/z direction (2D/3D)"
6363)
6464physical_parameters .add_argument (
65- "--gap" , type = float , default = - 0.05 , help = "y/z coordinate of rigid surface (2D/3D)"
65+ "--gap" , type = float , default = - 0.00 , help = "y/z coordinate of rigid surface (2D/3D)"
6666)
6767fem_parameters = parser .add_argument_group ("FEM parameters" )
6868fem_parameters .add_argument (
6969 "--degree" ,
7070 dest = "degree" ,
7171 type = int ,
72- default = 1 ,
72+ default = 2 ,
7373 help = "Degree of primal and latent space" ,
7474)
7575fem_parameters .add_argument (
@@ -81,7 +81,7 @@ class CustomParser(
8181 "--n-max-iterations" ,
8282 dest = "newton_max_iterations" ,
8383 type = int ,
84- default = 25 ,
84+ default = 250 ,
8585 help = "Maximum number of iterations of Newton iteration" ,
8686)
8787newton_parameters .add_argument (
@@ -113,7 +113,7 @@ class CustomParser(
113113alpha_options .add_argument (
114114 "--alpha_scheme" ,
115115 type = str ,
116- default = "constant " ,
116+ default = "doubling " ,
117117 choices = typing .get_args (AlphaScheme ),
118118 help = "Scheme for updating alpha" ,
119119)
@@ -240,7 +240,6 @@ def solve_contact_problem(
240240 # Define problem specific parameters
241241 mu = E / (2.0 * (1.0 + nu ))
242242 lmbda = E * nu / ((1.0 + nu ) * (1.0 - 2.0 * nu ))
243- n = ufl .FacetNormal (mesh )
244243 n_g = dolfinx .fem .Constant (mesh , np .zeros (gdim , dtype = dst ))
245244 n_g .value [- 1 ] = - 1
246245 alpha = dolfinx .fem .Constant (mesh , dst (alpha_0 ))
@@ -252,7 +251,7 @@ def solve_contact_problem(
252251 residual = alpha * ufl .inner (sigma (u , mu , lmbda ), epsilon (v )) * ufl .dx (
253252 domain = mesh
254253 ) - alpha * ufl .inner (f , v ) * ufl .dx (domain = mesh )
255- residual += - ufl .inner (psi - psi_k , ufl .dot (v , n )) * ds
254+ residual += - ufl .inner (psi - psi_k , ufl .dot (v , n_g )) * ds
256255 residual += ufl .inner (ufl .dot (u , n_g ), w ) * ds
257256 residual += ufl .inner (ufl .exp (psi ), w ) * ds - ufl .inner (g , w ) * ds
258257
@@ -371,7 +370,6 @@ def assemble_penetration():
371370 return it , iterations
372371
373372
374- # python3 script.py --alpha_0=0.1 --degree=2 file --filename=sphere.xdmf
375373if __name__ == "__main__" :
376374 args = parser .parse_args ()
377375 if args .mesh == "native" :
0 commit comments