Skip to content
Snippets Groups Projects
Commit 735749d3 authored by Sander, Oliver's avatar Sander, Oliver
Browse files

cosserat-continuum.cc: Move Dirichlet values into main problem file

That way, a single Python file is enough to describe the entire
boundary value problem.
parent 23329cca
No related branches found
No related tags found
1 merge request!159cosserat-continuum.cc: Read input parameters as Python file
Pipeline #15829 passed
import math
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
def deformation(self, x):
# Dirichlet b.c. simply clamp the shell in the reference configuration
out = [x[0], x[1], 0]
return out
def orientation(self, x):
rotation = [[1,0,0], [0, 1, 0], [0, 0, 1]]
return rotation
......@@ -106,13 +106,27 @@ parameterSet.materialParameters.b3 = 1
# Boundary values
#############################################
parameterSet.problem = "cantilever"
### Python predicate specifying all Dirichlet grid vertices
# x is the vertex coordinate
parameterSet.dirichletVerticesPredicate = "[x[0] < 0.01, x[0] < 0.01, x[0] < 0.01]"
parameterSet.dirichletRotationVerticesPredicate = "x[0] < 0.01"
### The actual Dirichlet values
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
def deformation(self, x):
# Dirichlet b.c. simply clamp the shell in the reference configuration
out = [x[0], x[1], 0]
return out
def orientation(self, x):
rotation = [[1,0,0], [0, 1, 0], [0, 0, 1]]
return rotation
### Python predicate specifying all Neumann grid vertices
# x is the vertex coordinate
parameterSet.neumannVerticesPredicate = "x[0] > 99.99"
......
......@@ -103,13 +103,38 @@ parameterSet.materialParameters.b3 = 1
# Boundary values
#############################################
parameterSet.problem = "twisted-strip"
### Python predicate specifying all Dirichlet grid vertices
# x is the vertex coordinate
parameterSet.dirichletVerticesPredicate = "[x[0] < 0.001 or x[0] > 0.0999, x[0] < 0.001 or x[0] > 0.0999, x[0] < 0.001 or x[0] > 0.0999]"
parameterSet.dirichletRotationVerticesPredicate = "x[0] < 0.001 or x[0] > 0.0999"
### The actual Dirichlet values
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
self.upper = [0.1, 0.01]
self.totalAngle = 6*math.pi
def deformation(self, x):
angle = self.totalAngle * x[0]/self.upper[0]
angle *= self.homotopyParameter
# Rotation matrix (around y-axis)
rotation = [[1,0,0], [0, math.cos(angle), -math.sin(angle)], [0, math.sin(angle), math.cos(angle)]]
# Matrix-vector product, vector is [x[0], x[1], 0]
out = [rotation[0][0]*x[0]+rotation[0][1]*x[1], rotation[1][0]*x[0]+rotation[1][1]*x[1], rotation[2][0]*x[0]+rotation[2][1]*x[1]]
return out
def orientation(self, x):
angle = self.totalAngle * x[0]/self.upper[0]
angle *= self.homotopyParameter
rotation = [[1,0,0], [0, math.cos(angle), -math.sin(angle)], [0, math.sin(angle), math.cos(angle)]]
return rotation
# Initial deformation
parameterSet.initialDeformation = "[x[0], x[1], 0]"
......
......@@ -102,13 +102,27 @@ parameterSet.materialParameters.b3 = 1
# Boundary values
#############################################
parameterSet.problem = "wong-pellegrino"
### Python predicate specifying all Dirichlet grid vertices
# x is the vertex coordinate
parameterSet.dirichletVerticesPredicate = "[x[1] < 0.0001 or x[1] > 0.128 - 0.0001, x[1] < 0.0001 or x[1] > 0.128 - 0.0001, x[1] < 0.0001 or x[1] > 0.128 - 0.0001]"
parameterSet.dirichletRotationVerticesPredicate = "x[1] < 0.0001 or x[1] > 0.128 - 0.0001"
### The actual Dirichlet values
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
def deformation(self, x):
out = [x[0], x[1], 0]
if x[1] > 0.128-1e-4 :
out[0] += 0.003 * self.homotopyParameter
out[1] += 0.0005
return out
def orientation(self, x):
rotation = [[1,0,0], [0, 1, 0], [0, 0, 1]]
return rotation
# Initial deformation
#parameterSet.startFromFile = True
parameterSet.initialIterateFilename = "cosserat_iterate_2.vtu"
......
......@@ -98,13 +98,31 @@ parameterSet.materialParameters.b3 = 1
# Boundary values
#############################################
parameterSet.problem = "wriggers-l-shape"
### Python predicate specifying all Dirichlet grid vertices
# x is the vertex coordinate
parameterSet.dirichletVerticesPredicate = "[x[0] < 1, x[0] < 1, x[0] < 1]"
parameterSet.dirichletRotationVerticesPredicate = "x[0] < 1"
### The actual Dirichlet values
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
# Deformation of 3d classical materials
def dirichletValues(self, x):
# Clamp the L-shape in its reference configuration
return [x[0], x[1], x[2]]
# Deformation of Cosserat shells
def deformation(self, x):
# Clamp the L-shape in its reference configuration
return [x[0], x[1], 0]
# Orientation of Cosserat materials
def orientation(self, x):
rotation = [[1,0,0], [0, 1, 0], [0, 0, 1]]
return rotation
### Python predicate specifying all Dirichlet grid vertices
# x is the vertex coordinate
parameterSet.neumannVerticesPredicate = "x[1] < -239"
......
import math
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
self.upper = [0.1, 0.01]
self.totalAngle = 6*math.pi
def deformation(self, x):
angle = self.totalAngle * x[0]/self.upper[0]
angle *= self.homotopyParameter
# Rotation matrix (around y-axis)
rotation = [[1,0,0], [0, math.cos(angle), -math.sin(angle)], [0, math.sin(angle), math.cos(angle)]]
# Matrix-vector product, vector is [x[0], x[1], 0]
out = [rotation[0][0]*x[0]+rotation[0][1]*x[1], rotation[1][0]*x[0]+rotation[1][1]*x[1], rotation[2][0]*x[0]+rotation[2][1]*x[1]]
return out
def orientation(self, x):
angle = self.totalAngle * x[0]/self.upper[0]
angle *= self.homotopyParameter
rotation = [[1,0,0], [0, math.cos(angle), -math.sin(angle)], [0, math.sin(angle), math.cos(angle)]]
return rotation
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
def deformation(self, x):
out = [x[0], x[1], 0]
if x[1] > 0.128-1e-4 :
out[0] += 0.003 * self.homotopyParameter
out[1] += 0.0005
return out
def orientation(self, x):
rotation = [[1,0,0], [0, 1, 0], [0, 0, 1]]
return rotation
class DirichletValues:
def __init__(self, homotopyParameter):
self.homotopyParameter = homotopyParameter
# Deformation of 3d classical materials
def dirichletValues(self, x):
# Clamp the L-shape in its reference configuration
return [x[0], x[1], x[2]]
# Deformation of Cosserat shells
def deformation(self, x):
# Clamp the L-shape in its reference configuration
return [x[0], x[1], 0]
# Orientation of Cosserat materials
def orientation(self, x):
rotation = [[1,0,0], [0, 1, 0], [0, 0, 1]]
return rotation
......@@ -471,12 +471,11 @@ int main (int argc, char *argv[]) try
// Set Dirichlet values
////////////////////////////////////////////////////////
Python::Reference dirichletValuesClass = Python::import(parameterSet.get<std::string>("problem") + "-dirichlet-values");
// Dirichlet boundary value function, depending on the homotopy parameter
Python::Callable dirichletValuesPythonClass = pyModule.get("DirichletValues");
Python::Callable C = dirichletValuesClass.get("DirichletValues");
// Call a constructor.
Python::Reference dirichletValuesPythonObject = C(homotopyParameter);
// Construct with a particular value of the homotopy parameter
Python::Reference dirichletValuesPythonObject = dirichletValuesPythonClass(homotopyParameter);
// Extract object member functions as Dune functions
auto deformationDirichletValues = Python::make_function<FieldVector<double,3> > (dirichletValuesPythonObject.get("deformation"));
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment