diff --git a/experiment/self-folding-box/self-folding-box_parameterBackup.py b/experiment/self-folding-box/self-folding-box_parameterBackup.py new file mode 100644 index 0000000000000000000000000000000000000000..2cb8a6d9cdf127ff97fe48712f7f9d038842e7c5 --- /dev/null +++ b/experiment/self-folding-box/self-folding-box_parameterBackup.py @@ -0,0 +1,542 @@ +import math +import numpy as np + +class ParameterSet(dict): + def __init__(self, *args, **kwargs): + super(ParameterSet, self).__init__(*args, **kwargs) + self.__dict__ = self + +parameterSet = ParameterSet() + + + +############################################# +# Paths +############################################# +parameterSet.resultPath = '/home/klaus/Desktop/Dune_bendIso/dune-microstructure/outputs_self-folding-box' +# parameterSet.outputPath = '/home/klaus/Desktop/Dune_bendIso/dune-microstructure/outputs_buckling_experiment' # This is currently still used in prestrainedMaterial +parameterSet.baseName= 'self-folding-box' + +############################################# +# Grid parameters +############################################# +# nX=8 +# nY=8 + +# parameterSet.structuredGrid = 'simplex' +# parameterSet.lower = '0 0' +# parameterSet.upper = '4 1' +# parameterSet.elements = str(nX)+' '+ str(nY) + +parameterSet.macroGridLevel = 4 #good results only for refinementLevel >=4 + + +parameterSet.structuredGrid = 'false' +parameterSet.gridPath = '/home/klaus/Desktop/Dune_bendIso/dune-microstructure/experiment/self-folding-box' +# parameterSet.gridFile = 'crossbox.msh' +# parameterSet.gridFile = 'crossbox_new.msh' +parameterSet.gridFile = 'Teest.msh' +# parameterSet.gridFile = 'example.msh' + +############################################# +# Options +############################################# +parameterSet.measure_analyticalError = False +parameterSet.measure_isometryError = False +parameterSet.vtkWrite_analyticalSurfaceNormal = False +parameterSet.vtkwrite_analyticalSolution = False + + +parameterSet.conforming_DiscreteJacobian = 0 +############################################# +# Solver parameters +############################################# +# Tolerance of the multigrid solver +parameterSet.tolerance = 1e-12 +# Maximum number of multigrid iterations +parameterSet.maxProximalNewtonSteps = 200 +# Initial regularization +parameterSet.initialRegularization = 200 +# Measure convergence +parameterSet.instrumented = 0 + + +#TEST RIemannian TR + +# parameterSet.Solver = "RiemannianTR" +# parameterSet.numIt = 200 +# parameterSet.nu1 = 3 +# # Number of postsmoothing steps +# parameterSet.nu2 = 3 +# # Number of coarse grid corrections +# parameterSet.mu = 1 +# # Number of base solver iterations +# parameterSet.baseIt = 100 +# # Tolerance of the multigrid solver +# parameterSet.mgTolerance = 1e-10 +# # Tolerance of the base grid solver +# parameterSet.baseTolerance = 1e-8 +# parameterSet.tolerance = 1e-12 +# # Max number of steps of the trust region solver +# parameterSet.maxTrustRegionSteps = 100 +# # Initial trust-region radius +# parameterSet.initialTrustRegionRadius = 1 + + +############################ +# Problem specifications +############################ +# Dimension of the domain (only used for numerical-test python files) +parameterSet.dim = 2 + +############################################# +# VTK/Output +############################################# +# Write discrete solution as .vtk-File +parameterSet.writeVTK = 1 +parameterSet.vtkwrite_analyticalSolution = 0 +parameterSet.vtkWrite_analyticalSurfaceNormal = 0 + +# The grid can be refined several times for a higher resolution in the VTK-file. +parameterSet.subsamplingRefinement = 0 + +# Write Dof-Vector to .txt-file +parameterSet.writeDOFvector = 0 + +############################################# +# Dirichlet boundary indicator +############################################# +# def one_norm(v,w): +# return np.abs(v) + np.abs(w) + + +def dirichlet_indicator(x) : + midpoint = [0.5,2.5] + # if(one_norm(x[0]-midpoint[0],x[1]-midpoint[1]) <= 0.5 ): + if(np.max([abs(x[0]-midpoint[0]),abs(x[1]-midpoint[1])]) <= 0.5 ): + return True + else: + return False + + +# def dirichlet_indicator(x) : +# # if(((x[0] > -0.01) and (x[0] < 1.01)) and ((x[1] > 1.99) and (x[1] < 3.01))): +# return True +# else: +# return False + +#TEST +# def dirichlet_indicator(x): +# if( x[1] <= 1.0): +# return True +# else: +# return False + + + + +############################################# +# Initial iterate function +############################################# +def f(x): + return [x[0], x[1], 0] + +def df(x): + return ((1,0), + (0,1), + (0,0)) +# #Rotation: +# def R(beta): +# return [[math.cos(beta),0], +# [0,1], +# [-math.sin(beta),0]] + + +# def f(x): +# a = 0.5 +# if(x[0] <= 0.01): +# return [x[0], x[1], 0] +# elif(x[0] >= 3.99): +# return [x[0] - a, x[1], 0] +# else: +# return [x[0], x[1], 0] + + +# # beta = math.pi/4.0 +# # beta= 0.05 +# beta = 0 + +# def df(x): +# a = 0.5 +# if(x[0] <= 0.01): +# # return R(-1.0*beta) +# return R(beta) +# elif(x[0] >= 3.99): +# # return R(beta) +# return R(-1.0*beta) +# else: +# return ((1,0), +# (0,1), +# (0,0)) + + + + +fdf = (f, df) + + +############################################# +# Force +############################################ +parameterSet.assemble_force_term = False + +def force(x): + return [0, 0, 0] + + + + +##################### MICROSCALE PROBLEM #################### + +# parameterSet.prestrainFlag = True #deprecated +# parameterSet.macroscopically_varying_microstructure = False +# parameterSet.read_effectiveQuantities_from_Parset = True # Otherwise the Micro/Cell-Problem is solved once to obtain the quantities. +parameterSet.printMicroOutput = False + +#New +# parameterSet.piecewiseConstantMicrostructure = True +parameterSet.VTKwriteMacroPhaseIndicator = True +parameterSet.MacroPhaseSubsamplingRefinement = 3 + + + + + +class GlobalMicrostructure: + """ Class that represents the global microstructure. + + The global microstructure can be defined individually + for different (finitely many) macroscopic phases. + Each macroscopic phase corresponds to a 'LocalMicrostructure_' + sub-class that defines the necessary data for the local + micro-problem. + + Currently, there are three possibilities for the LocalMicrostructure: + (1) Read the effective quantities (Qhom,Beff) from this parset. + (Constant local microstructure) + (2) Solve the micro-problem once to obtain the effective quantities. + (Constant local microstructure) + (3) Solve for micro-problem for each macroscopic quadrature point. + (Macroscopocally varying local microstructure)) + """ + def __init__(self,macroPoint=[0,0]): + self.macroPhases = 4 + + # define first local microstructure options. + self.macroPhase1_constantMicrostructure = True + self.macroPhase1_readEffectiveQuantities = False; + + # define second local microstructure options. + self.macroPhase2_constantMicrostructure = True + self.macroPhase2_readEffectiveQuantities = False; + + # define third local microstructure options. + self.macroPhase3_constantMicrostructure = True + self.macroPhase3_readEffectiveQuantities = True; + + + self.macroPhase4_constantMicrostructure = True + self.macroPhase4_readEffectiveQuantities = True; + + + # def macroPhaseIndicator(self,y): #y :macroscopic point + # """ Indicatorfunction that determines the domains + # i.e. macro-phases of different local microstructures. + # """ + # a = 1.0/8.0 + # # self.macroPoint = x #just for visualization purposes + # if(y[1] < 1.0 +a and y[1] > 1.0): + # return 1 + # elif(y[1] < 2.0 + a and y[1] > 2.0): + # return 1 + # elif(y[1] < 3.0 + a and y[1] > 3.0): + # return 1 + # elif(y[0] < 0.0 and y[0] > -a): + # return 2 + # elif(y[0] < 1.0 + a and y[0] > 1.0): + # return 2 + # else: + # return 3 + + #old version: middle Cell clamped. + def macroPhaseIndicator(self,y): #y :macroscopic point + """ Indicatorfunction that determines the domains + i.e. macro-phases of different local microstructures. + """ + a = 1.0/20.0 + # self.macroPoint = x #just for visualization purposes + if(y[1] < 1.0 and y[1] > 1.0-a): + return 4 + elif(y[1] < 2.0 and y[1] > 2.0-a): + return 1 + elif(y[1] < 3.0 + a and y[1] > 3.0): + return 1 + elif(y[0] < 0.0 and y[0] > -a): + return 2 + elif(y[0] < 1.0 + a and y[0] > 1.0): + return 2 + else: + return 3 + + # if(y[1] < 2.0 and y[1] > 2.0-a): + # return 1 + # elif(y[1] < 3.0 + a and y[1] > 3.0): + # return 1 + # elif(y[0] < 0.0 and y[0] > -a): + # return 2 + # elif(y[0] < 1.0 + a and y[0] > 1.0): + # return 2 + # else: + # return 3 + + + + # Represents the local microstructure in Phase 1 + class LocalMicrostructure_1: + def __init__(self,macroPoint=[0,0]): + self.macroPoint = macroPoint + self.gamma = 1.0 #in the future this might change depending on macroPoint. + self.phases = 2 #in the future this might change depending on macroPoint. + #--- Define different local material phases: + self.phase1_type="isotropic" + self.materialParameters_phase1 = [200, 1.0] + + self.phase2_type="isotropic" + self.materialParameters_phase2 = [100, 1.0] + + + #--- Indicator function for material phases + def indicatorFunction(self,x): + # if (abs(x[0]) < (theta/2) and x[2] >= 0 ): + # if (abs(x[0]) < (1.0/4.0) and x[2] <= 0 ): + if (x[2] <= 0 ): + return 1 #Phase1 + else : + return 2 #Phase2 + + #--- Define prestrain function for each phase (also works with non-constant values) + def prestrain_phase1(self,x): + # return [[2, 0, 0], [0,2,0], [0,0,2]] + # rho = 0 + rho = 25.0 + return [[rho*1.0, 0, 0], [0,rho*1.0,0], [0,0,rho*1.0]] + + def prestrain_phase2(self,x): + return [[0, 0, 0], [0,0,0], [0,0,0]] + + + # Represents the local microstructure in Phase 1 + class LocalMicrostructure_2: + def __init__(self,macroPoint=[0,0]): + self.macroPoint = macroPoint + self.gamma = 1.0 #in the future this might change depending on macroPoint. + self.phases = 2 #in the future this might change depending on macroPoint. + #--- Define different local material phases: + self.phase1_type="isotropic" + self.materialParameters_phase1 = [200, 1.0] + + self.phase2_type="isotropic" + self.materialParameters_phase2 = [100, 1.0] + + + #--- Indicator function for material phases + def indicatorFunction(self,x): + # if (abs(x[0]) < (theta/2) and x[2] >= 0 ): + # if (abs(x[1]) < (1.0/4.0) and x[2] <= 0 ): + if (x[2] <= 0 ): + return 1 #Phase1 + else : + return 2 #Phase2 + + #--- Define prestrain function for each phase (also works with non-constant values) + def prestrain_phase1(self,x): + # return [[2, 0, 0], [0,2,0], [0,0,2]] + # rho = 0 + rho = 20.0 + return [[rho*1.0, 0, 0], [0,rho*1.0,0], [0,0,rho*1.0]] + + def prestrain_phase2(self,x): + return [[0, 0, 0], [0,0,0], [0,0,0]] + + + # Represents the local microstructure in Phase 2 + class LocalMicrostructure_3: + def __init__(self,macroPoint=[0,0]): + self.macroPoint = macroPoint + self.gamma = 1.0 + + self.effectivePrestrain= np.array([[0.0, 0.0], + [0.0, 0.0]]) + + self.effectiveQuadraticForm = np.array([[1.0, 0.0, 0.0], + [0.0, 1.0, 0.0], + [0.0, 0.0, 1.0]]) + + + # Represents the local microstructure in Phase 1 + class LocalMicrostructure_4: + def __init__(self,macroPoint=[0,0]): + self.macroPoint = macroPoint + self.gamma = 1.0 #in the future this might change depending on macroPoint. + self.phases = 2 #in the future this might change depending on macroPoint. + #--- Define different local material phases: + self.phase1_type="isotropic" + self.materialParameters_phase1 = [200, 1.0] + + self.phase2_type="isotropic" + self.materialParameters_phase2 = [100, 1.0] + + + self.effectivePrestrain= np.array([[0.0, 0.0], + [0.0, -10.0]]) + + self.effectiveQuadraticForm = np.array([[23.0, 0.0898, 0.0], + [0.0898, 23.0, 0.0], + [0.0, 0.0, 22.9]]) + + + + #--- Indicator function for material phases + def indicatorFunction(self,x): + # if (abs(x[0]) < (theta/2) and x[2] >= 0 ): + # if (abs(x[0]) < (1.0/4.0) and x[2] <= 0 ): + # if (x[2] <= 0 ): + if (x[1] <= 0 ): + return 1 #Phase1 + else : + return 2 #Phase2 + + #--- Define prestrain function for each phase (also works with non-constant values) + def prestrain_phase1(self,x): + # return [[2, 0, 0], [0,2,0], [0,0,2]] + # rho = 0 + rho = 8.0 + return [[rho*1.0, 0, 0], [0,rho*1.0,0], [0,0,rho*1.0]] + + def prestrain_phase2(self,x): + return [[0, 0, 0], [0,0,0], [0,0,0]] + + + +# # Represents the local microstructure +# class Microstructure: +# def __init__(self,macroPoint=[0,0]): #default value for initialization +# self.macroPoint = macroPoint +# # self.macroPoint = macroPoint +# self.gamma = 1.0 #in the future this might change depending on macroPoint. +# self.phases = 2 #in the future this might change depending on macroPoint. +# #--- Define different material phases: +# #- PHASE 1 +# self.phase1_type="isotropic" +# self.materialParameters_phase1 = [200, 1.0] +# #- PHASE 2 +# self.phase2_type="isotropic" +# self.materialParameters_phase2 = [100, 1.0] + +# self.effectivePrestrain= np.array([[1.0, 0.0], +# [0.0, 1.0]]) + +# self.effectiveQuadraticForm = np.array([[1.0, 0.0, 0.0], +# [0.0, 1.0, 0.0], +# [0.0, 0.0, 1.0]]) + +# # --------------------------- +# self.macroPhases = 2 + +# self.macroPhase1_readEffectiveQuantities = False; + + + + +# def macroPhaseIndicator(self,x): +# a = 1.0/4.0 + +# self.macroPoint = x #just for visualization purposes + +# if(self.macroPoint[1] < 1.0 and self.macroPoint[1] > 1.0-a): +# return 1 +# elif(self.macroPoint[1] < 2.0 and self.macroPoint[1] > 2.0-a): +# return 1 +# elif(self.macroPoint[1] < 3.0 + a and self.macroPoint[1] > 3.0): +# return 1 +# elif(self.macroPoint[0] < 0.0 and self.macroPoint[0] > -a): +# return 1 +# elif(self.macroPoint[0] < 1.0 + a and self.macroPoint[0] > 1.0): +# return 1 +# else: +# return 0 + + +# #--- Indicator function for material phases +# def indicatorFunction(self,x): +# # if (abs(x[0]) < (theta/2) and x[2] >= 0 ): +# if (abs(x[0]) < (1.0/4.0) and x[2] <= 0 ): +# return 1 #Phase1 +# else : +# return 2 #Phase2 + +# #--- Define prestrain function for each phase (also works with non-constant values) +# def prestrain_phase1(self,x): +# # return [[2, 0, 0], [0,2,0], [0,0,2]] +# # rho = 0 +# rho = 1.0 +# return [[rho*1.0, 0, 0], [0,rho*1.0,0], [0,0,rho*1.0]] + +# def prestrain_phase2(self,x): +# return [[0, 0, 0], [0,0,0], [0,0,0]] + +# + + + +############################################# +# Grid parameters +############################################# +parameterSet.microGridLevel = 3 + +############################################# +# Assembly options +############################################# +parameterSet.set_IntegralZero = 1 #(default = false) +parameterSet.set_oneBasisFunction_Zero = 1 #(default = false) +#parameterSet.arbitraryLocalIndex = 7 #(default = 0) +#parameterSet.arbitraryElementNumber = 3 #(default = 0) + +############################################# +# Solver Options, Type: #1: CG - SOLVER , #2: GMRES - SOLVER, #3: QR - SOLVER, #4: UMFPACK - SOLVER (default) +############################################# +parameterSet.Solvertype = 4 # recommended to use iterative solver (e.g GMRES) for finer grid-levels +parameterSet.Solver_verbosity = 0 #(default = 2) degree of information for solver output + + +############################################# +# Write/Output options #(default=false) +############################################# +# --- (Optional output) write Material / prestrain / Corrector functions to .vtk-Files: +parameterSet.write_materialFunctions = 1 # VTK indicator function for material/prestrain definition +#parameterSet.write_prestrainFunctions = 1 # VTK norm of B (currently not implemented) +parameterSet.MaterialSubsamplingRefinement= 2 + +# --- (Additional debug output) +parameterSet.print_debug = 0 #(default=false) + +parameterSet.print_corrector_matrices = 0 + +# --- Write Correctos to VTK-File: +parameterSet.writeCorrectorsVTK = 1 + +# --- Use caching of element matrices: +parameterSet.cacheElementMatrices = 1 + +# --- check orthogonality (75) from paper: +parameterSet.write_checkOrthogonality = 0