strain softening

docker/mandyoc/Dockerfile

@@ -32,6 +32,7 @@ COPY --chown=${USER}:${USER} . ${MANDYOC_DIR}
 RUN mkdir -p ${HOME}/.local/bin && \
     cd ${MANDYOC_DIR} && \
     make clear && \
+    make all && \
     make install
 
 # =============================================================================

src/DMSwarm_2d.cpp

@@ -64,6 +64,14 @@ extern PetscInt seed_layer_size;
 extern PetscBool seed_layer_set;
 extern PetscBool strain_seed_constant;
 
+extern PetscBool weakening_from_interfaces_file;
+
+// From param.txt and interfaces.txt
+extern PetscInt WITH_NON_LINEAR;
+extern PetscInt PLASTICITY;
+extern PetscScalar *weakening_seed;
+extern PetscScalar weakening_min;
+
 extern PetscReal *strain_seed_layer;
 extern PetscInt strain_seed_layer_size;
 extern PetscBool strain_seed_layer_set;
@@ -442,31 +450,15 @@ PetscErrorCode createSwarm_2d()
 					layer_array[p] = n_interfaces;
 				}
 
+				// Load initial strain_array with weakening_seed to or random value
 				rand_r(&seed_strain);
-				strain_array[p]=random_initial_strain*(float)rand_r(&seed_strain)/RAND_MAX;
-
-				if (seed_layer_set == PETSC_TRUE) {
-					if (seed_layer_size == 1 && layer_array[p] == seed_layer[0]) {
-						if (strain_seed_constant == PETSC_TRUE) {
-							strain_array[p] = strain_seed_layer[0];
-						} else {
-							strain_array[p] += strain_seed_layer[0];
-						}
-					}
-					else {
-						for (int k = 0; k < seed_layer_size; k++) {
-							if (layer_array[p] == seed_layer[k]) {
-								if (strain_seed_constant == PETSC_TRUE) {
-									strain_array[p] = strain_seed_layer[k];
-								} else {
-									strain_array[p] += strain_seed_layer[k];
-								}
-							}
-						}
-					}
+				if (weakening_seed[layer_array[p]] >= 0) {
+					strain_array[p] = weakening_seed[layer_array[p]];
+				} 
+				else
+				{
+					strain_array[p] = random_initial_strain*(float)rand_r(&seed_strain)/RAND_MAX;
 				}
-
-
 			}
 		}
 

src/DMSwarm_3d.cpp

@@ -65,6 +65,10 @@ extern PetscInt print_step_files;
 
 extern PetscReal random_initial_strain;
 
+// From param.txt and interfaces.txt
+extern PetscInt WITH_NON_LINEAR;
+extern PetscInt PLASTICITY;
+extern PetscScalar *weakening_seed;
 
 extern PetscInt *seed_layer;
 extern PetscInt seed_layer_size;
@@ -475,28 +479,14 @@ PetscErrorCode createSwarm_3d()
 					//printf("entrei!\n");
 				}
 
+				// Load initial strain_array with weakening_seed to or random value
 				rand_r(&seed_strain);
-				strain_array[p]=random_initial_strain*(float)rand_r(&seed_strain)/RAND_MAX;
-
-				if (seed_layer_set == PETSC_TRUE) {
-					if (seed_layer_size == 1 && layer_array[p] == seed_layer[0]) {
-						if (strain_seed_constant == PETSC_TRUE) {
-							strain_array[p] = strain_seed_layer[0];
-						} else {
-							strain_array[p] += strain_seed_layer[0];
-						}
-					}
-					else {
-						for (int k = 0; k < seed_layer_size; k++) {
-							if (layer_array[p] == seed_layer[k]) {
-								if (strain_seed_constant == PETSC_TRUE) {
-									strain_array[p] = strain_seed_layer[k];
-								} else {
-									strain_array[p] += strain_seed_layer[k];
-								}
-							}
-						}
-					}
+				if (weakening_seed[layer_array[p]] >= 0) {
+					strain_array[p] = weakening_seed[layer_array[p]];
+				} 
+				else
+				{
+					strain_array[p] = random_initial_strain*(float)rand_r(&seed_strain)/RAND_MAX;
 				}
 				/////!!!!
 				/*if (rank==0){

src/DMSwarm_move.cpp

@@ -18,7 +18,7 @@ typedef struct {
 } Stokes3d;
 
 double calc_visco_ponto(double T, double P, double x, double z,double geoq_ponto,double e2_inva,double strain_cumulate,
-						double A, double n_exp, double QE, double VE);
+						double A, double n_exp, double QE, double VE, PetscInt layer_number);
 
 extern DM dms;
 
@@ -383,7 +383,8 @@ PetscErrorCode moveSwarm(int dimensions, PetscReal dt)
 
 
 			rarray[p] = calc_visco_ponto(tp,Pp,cx,cz,inter_geoq[layer_array[p]],E2_invariant,strain_fac[p],
-										inter_A[layer_array[p]], inter_n[layer_array[p]], inter_Q[layer_array[p]], inter_V[layer_array[p]]);
+										inter_A[layer_array[p]], inter_n[layer_array[p]], inter_Q[layer_array[p]], inter_V[layer_array[p]], layer_array[p]);
+			// PetscPrintf(PETSC_COMM_WORLD, "p: %d, strain: %E\n", p, strain_fac[p]);
 
 
 
@@ -575,8 +576,10 @@ PetscErrorCode moveSwarm(int dimensions, PetscReal dt)
 
 
 			rarray[p] = calc_visco_ponto(tp,Pp,cx,cz,inter_geoq[layer_array[p]],E2_invariant,strain_fac[p],
-										inter_A[layer_array[p]], inter_n[layer_array[p]], inter_Q[layer_array[p]], inter_V[layer_array[p]]);
+										inter_A[layer_array[p]], inter_n[layer_array[p]], inter_Q[layer_array[p]], inter_V[layer_array[p]], layer_array[p]);
 
+			// PetscPrintf(PETSC_COMM_WORLD, "p: %d, strain: %E\n", p, strain_fac[p]);
+
 			array[3*p  ] += dt * vx;
 			array[3*p+1] += dt * vy;
 			array[3*p+2] += dt * vz;