rename geometries-folder to materials

materials/material_backup.py

+import math
+from python_matrix_operations import *
+import ctypes
+import os
+import sys
+
+
+Phases=3
+
+mu_ = [80, 80, 60]
+lambda_ = [80, 80, 25]
+
+
+phase1_type="isotropic"
+materialParameters_phase1 = [80, 80]
+
+phase2_type="isotropic"
+materialParameters_phase2 = [80, 80]
+
+phase3_type="isotropic"
+materialParameters_phase3 = [60, 25]
+
+
+# print('mu_:', mu_)
+# A = [[1, 5, 0], [5,1,0], [5,0,1]]
+#
+# print("sym(A):", sym(A))
+
+
+# dir_path = os.path.dirname(os.path.realpath("/home/klaus/Desktop/Dune-Testing/dune-microstructure/dune/microstructure/matrix_operations.hh"))
+# handle = ctypes.CDLL(dir_path)
+#
+# handle.create2Darray.argtypes = [ctypes.c_int, ctypes.c_double, ctypes.c_double]
+#
+# def create2Darray(nside, mx, my):
+#     return handle.create2Darray(nside, mx, my)
+
+
+#Indicator function that determines both phases
+# x[0] : y1-component
+# x[1] : y2-component
+# x[2] : x3-component
+#    --- replace with your definition of indicatorFunction:
+###############
+# Wood
+###############
+# def f(x):
+#     theta=0.25
+#     # mu_ = [3, 5, 6]
+#     # lambda_ = [9, 7, 8]
+#     # mu_ = 3 5 6
+#     # lambda_ = 9 7 8
+
+#     if ((abs(x[0]) < theta/2) and x[2]<0.25):
+#         return [mu_[0], lambda_[0]]    #latewood
+#         # return 5    #latewood
+#     elif ((abs(x[0]) > theta/2) and x[2]<0.25):
+#         return [mu_[1], lambda_[1]]    #latewood
+#         # return 2
+#     else :
+#         return [mu_[2],lambda_[2]]    #latewood    #Phase3
+#         # return 1
+
+def indicatorFunction(x):
+    theta=0.25
+    factor=1
+    if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+        return 1    #Phase1
+    elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+        return 2    #Phase2
+    else :
+        return 0    #Phase3
+
+
+def L1(G):
+    return 2.0 * mu_[0] * sym(G) + lambda_[0] * trace(sym(G)) * Id()   #Phase1
+
+def L2(G):
+    return 2.0 * mu_[0] * sym(G) + lambda_[0] * trace(sym(G)) * Id()   #Phase1
+
+def L3(G):
+    return 2.0 * mu_[0] * sym(G) + lambda_[0] * trace(sym(G)) * Id()   #Phase1
+
+
+# TEST
+
+# def L1(G):
+#     return Add(smult(2.0 * mu_[0], sym(G)),smult(lambda_[0] ,smult(trace(sym(G)),Id()) ))  #Phase1
+
+# def L2(G):
+#     return Add(smult(2.0 * mu_[1], sym(G)),smult(lambda_[1] ,smult(trace(sym(G)),Id()) ))   #Phase1
+
+# def L3(G):
+#     return Add(smult(2.0 * mu_[2], sym(G)),smult(lambda_[2] ,smult(trace(sym(G)),Id()) ))   #Phase1
+
+
+#Workaround
+def muValue(x):
+    return mu_
+
+def lambdaValue(x):
+    return lambda_
+
+
+
+
+
+###############
+# Cross
+###############
+# def f(x):
+#     theta=0.25
+#     factor=1
+#     if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+#         return 1    #Phase1
+#     elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+#         return 2    #Phase2
+#     else :
+#         return 0    #Phase3
+
+
+
+
+
+# def f(x):
+#     # --- replace with your definition of indicatorFunction:
+#     if (abs(x[0]) > 0.25):
+#         return 1    #Phase1
+#     else :
+#         return 0    #Phase2
+
+def b1(x):
+    return [[1, 0, 0], [0,1,0], [0,0,1]]
+
+def b2(x):
+    return [[1, 0, 0], [0,1,0], [0,0,1]]
+
+def b3(x):
+    return [[0, 0, 0], [0,0,0], [0,0,0]]
+
+# mu=80 80 60
+
+# lambda=80 80 25
+
+
+# --- elasticity tensor
+# def L(G,x):
+# def L(G):
+#     # input:  symmetric matrix G, position x
+#     # output: symmetric matrix LG
+#     return [[1, 1, 1], [1, 1, 1],[1, 1, 1]]
+
+
+
+
+
+
+
+
+
+# --- elasticity tensor
+def L(G,x):
+    # input:  symmetric matrix G, position x
+    # output: symmetric matrix LG
+    theta=0.25
+    if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+        return 2.0 * mu_[0] * sym(G) + lambda_[0] * trace(sym(G)) * Id()   #Phase1
+    elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+        return 2.0 * mu_[1] * sym(G) + lambda_[1] * trace(sym(G)) * Id()   #Phase2
+    else :
+        return 2.0 * mu_[2] * sym(G) + lambda_[2] * trace(sym(G)) * Id()   #Phase3
+# #     # 2.0 * mu * sym(E1) + lambda * trace(sym(E1)) * Id();
+
+
+# def L(G,x):
+#     # input:  symmetric matrix G, position x
+#     # output: symmetric matrix LG
+#     theta=0.25
+#     if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+#         # return 2.0 * mu_[0] * sym(G) + lambda_[0] * trace(sym(G)) * Id()   #Phase1
+#         return Add(smult(2.0 * mu_[0], sym(G)),smult(lambda_[0] ,smult(trace(sym(G)),Id()) ))
+#     elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+#         return Add(smult(2.0 * mu_[0], sym(G)),smult(lambda_[0] ,smult(trace(sym(G)),Id()) ))   #Phase2
+#     else :
+#         return Add(smult(2.0 * mu_[0], sym(G)),smult(lambda_[0] ,smult(trace(sym(G)),Id()) ))   #Phase3
+#         # return [[0, 0, 0], [0,0,0], [0,0,0]]  #Phase3
+
+
+##TEST
+# def L(G,x):
+#     # input:  symmetric matrix G, position x
+#     # output: symmetric matrix LG
+#     theta=0.25
+#     if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+#         return [[1, 1, 1], [1, 1, 1],[1, 1, 1]]
+#     elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+#         return [[x[0], 1, x[0]], [1, 1, 1],[x[0], x[0], x[0]]]
+#     else :
+#         return [[0, x[2], x[2]], [0,x[2],0], [0,0,0]]
+
+
+##TEST
+# def L(G,x):
+#     # input:  symmetric matrix G, position x
+#     # output: symmetric matrix LG
+#     theta=0.25
+#     if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+#         return sym([[1, 1, 1], [1, 1, 1],[1, 1, 1]])
+#     elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+#         return sym([[x[0], 1, x[0]], [1, 1, 1],[x[0], x[0], x[0]]])
+#     else :
+#         return sym([[0, x[2], x[2]], [0,x[2],0], [0,0,0]])
+
+
+
+# # small speedup..
+# def L(G,x):
+#     # input:  symmetric matrix G, position x
+#     # output: symmetric matrix LG
+#     theta=0.25
+#     if (x[0] <-1/2+theta and x[2]<-1/2+theta):
+#         return mu_[0] * (np.array(G).transpose() + np.array(G)) + lambda_[0] * (G[0][0] + G[1][1] + G[2][2]) * np.identity(3)   #Phase1
+#     elif (x[1]< -1/2+theta and x[2]>1/2-theta):
+#         return mu_[1] * (np.array(G).transpose() + np.array(G)) + lambda_[1] * (G[0][0] + G[1][1] + G[2][2]) * np.identity(3)   #Phase2
+#     else :
+#         return mu_[2] * (np.array(G).transpose() + np.array(G)) + lambda_[2] * (G[0][0] + G[1][1] + G[2][2]) * np.identity(3)   #Phase3
+#     # 2.0 * mu * sym(E1) + lambda * trace(sym(E1)) * Id();
+
+
+
+
+
+
+# def H(G,x):
+#     # input:  symmetric matrix G, position x
+#     # output: symmetric matrix LG
+#     if (abs(x[0]) > 0.25):
+#         return [[1, 1, 1], [1, 1, 1],[1, 1, 1]]
+#     else:
+#         return [[0, 0, 0], [0,0,0], [0,0,0]]
+
+
+def H(G,x):
+    # input:  symmetric matrix G, position x
+    # output: symmetric matrix LG
+    if (abs(x[0]) > 0.25):
+        return [[1, 1, 1], [1, 1, 1],[1, 1, 1]]
+    else:
+        return [[0, 0, 0], [0,0,0], [0,0,0]]
+
+# 2.0 * mu * sym(E1) + lambda * trace(sym(E1)) * Id();

geometries/material.py → materials/material_test.py

@@ -50,7 +50,6 @@ materialParameters_phase1 = [80, 80]
 phase2_type="isotropic"
 materialParameters_phase2 = [80, 80]
 
-
 #- PHASE 3
 phase3_type="isotropic"
 materialParameters_phase3 = [60, 25]
@@ -65,7 +64,8 @@ materialParameters_phase3 = [60, 25]
 #                                       [0.0,     0.0,     0.0,   0.0 ,         0.0,   1.0]])
 
 
-# define prestrain function for each phase
+#--- define prestrain function for each phase
+# (also works with non-constant values)
 def prestrain_phase1(x):
     return [[1, 0, 0], [0,1,0], [0,0,1]]