Commit c89c51e4 authored by Baptiste Durand's avatar Baptiste Durand

Remove periodic boundary condition for strain function space

Periodic BCs do not seem relevant for DG function spaces. In addition they may not be supported for DG function spaces for now : https://fenicsproject.discourse.group/t/dg-for-advection-equation-with-periodic-boundary-condition/272
Test for demo_homog2D.py : values of Chom and D# have not changed.
+ Black
parent 4ad7f8b0
No preview for this file type
......@@ -6,23 +6,27 @@ import numpy as np
from ho_homog import periodicity
from ho_homog.materials import epsilon, sigma, cross_energy
from ho_homog.toolbox_FEniCS import (DOLFIN_KRYLOV_METHODS, DOLFIN_LU_METHODS,
local_project)
from ho_homog.toolbox_FEniCS import (
DOLFIN_KRYLOV_METHODS,
DOLFIN_LU_METHODS,
local_project,
)
SOLVER_METHOD = 'mumps'
SOLVER_METHOD = "mumps"
np.set_printoptions(precision=4, linewidth=150)
np.set_printoptions(suppress=True)
logging.getLogger('UFL').setLevel(logging.DEBUG)
logging.getLogger('FFC').setLevel(logging.DEBUG)
logging.getLogger("UFL").setLevel(logging.DEBUG)
logging.getLogger("FFC").setLevel(logging.DEBUG)
# https://fenicsproject.org/qa/3669/how-can-i-diable-log-and-debug-messages-from-ffc/
logger = logging.getLogger(__name__)
class Fenics2DHomogenization(object):
''' Homogenization of 2D orthotropic periodic laminates in plane strain. The direction of lamination is 3 and the invariant direction is 1
'''
""" Homogenization of 2D orthotropic periodic laminates in plane strain. The direction of lamination is 3 and the invariant direction is 1
"""
def __init__(self, fenics_2d_rve, **kwargs):
"""[summary]
......@@ -45,13 +49,14 @@ class Fenics2DHomogenization(object):
self.rve = fenics_2d_rve
self.topo_dim = topo_dim = fenics_2d_rve.mesh_dim
self.pbc = periodicity.PeriodicDomain.pbc_dual_base(
fenics_2d_rve.gen_vect, 'XY', topo_dim)
fenics_2d_rve.gen_vect, "XY", topo_dim
)
solver = kwargs.pop('solver', {})
solver = kwargs.pop("solver", {})
# {'type': solver_type, 'method': solver_method, 'preconditioner': preconditioner}
s_type = solver.pop('type', None)
s_method = solver.pop('method', SOLVER_METHOD)
s_precond = solver.pop('preconditioner', None)
s_type = solver.pop("type", None)
s_method = solver.pop("method", SOLVER_METHOD)
s_precond = solver.pop("preconditioner", None)
if s_type is None:
if s_method in DOLFIN_KRYLOV_METHODS.keys():
s_type = "Krylov"
......@@ -61,11 +66,11 @@ class Fenics2DHomogenization(object):
raise RuntimeError("The indicated solver method is unknown.")
self._solver = dict(type=s_type, method=s_method)
if s_precond:
self._solver['preconditioner'] = s_precond
self._solver["preconditioner"] = s_precond
element = kwargs.pop('element', ('Lagrange', 2))
element = kwargs.pop("element", ("Lagrange", 2))
if isinstance(element, dict):
element = (element['family'], element['degree'])
element = (element["family"], element["degree"])
self._element = element
# * Function spaces
......@@ -74,29 +79,29 @@ class Fenics2DHomogenization(object):
self.displ_FE = fe.VectorElement(element[0], cell, element[1])
strain_deg = element[1] - 1 if element[1] >= 1 else 0
strain_dim = int(topo_dim * (topo_dim + 1) / 2)
self.strain_FE = fe.VectorElement('DG', cell, strain_deg,
dim=strain_dim)
self.strain_FE = fe.VectorElement("DG", cell, strain_deg, dim=strain_dim)
# Espace fonctionel scalaire
self.X = fe.FunctionSpace(self.rve.mesh, self.scalar_FE,
constrained_domain=self.pbc)
self.X = fe.FunctionSpace(
self.rve.mesh, self.scalar_FE, constrained_domain=self.pbc
)
# Espace fonctionnel 3D : deformations, notations de Voigt
self.W = fe.FunctionSpace(self.rve.mesh, self.strain_FE,
constrained_domain=self.pbc)
self.W = fe.FunctionSpace(self.rve.mesh, self.strain_FE)
# Espace fonctionel 2D pour les champs de deplacement
# TODO : reprendre le Ve défini pour l'espace fonctionnel mixte. Par ex: V = FunctionSpace(mesh, Ve)
self.V = fe.VectorFunctionSpace(self.rve.mesh, element[0], element[1],
constrained_domain=self.pbc)
self.V = fe.VectorFunctionSpace(
self.rve.mesh, element[0], element[1], constrained_domain=self.pbc
)
# * Espace fonctionel mixte pour la résolution : 2D pour les champs + scalaire pour multiplicateur de Lagrange
# * Espace fonctionel mixte pour la résolution :
# * 2D pour les champs + scalaire pour multiplicateur de Lagrange
# Pour le multiplicateur de Lagrange : Real element with one global degree of freedom
# "R" : Real element with one global degree of freedom
self.real_FE = fe.VectorElement("R", cell, 0)
self.M = fe.FunctionSpace(
self.rve.mesh,
fe.MixedElement([self.displ_FE, self.real_FE]),
constrained_domain=self.pbc
)
constrained_domain=self.pbc,
)
# Define variational problem
self.v, self.lamb_ = fe.TestFunctions(self.M)
......@@ -104,24 +109,28 @@ class Fenics2DHomogenization(object):
self.w = fe.Function(self.M)
# bilinear form
self.a = (fe.inner(
sigma(self.rve.C_per, epsilon(self.u)),
epsilon(self.v)
) * fe.dx
+ fe.dot(self.lamb_, self.u)*fe.dx
+ fe.dot(self.lamb, self.v)*fe.dx)
self.a = (
fe.inner(sigma(self.rve.C_per, epsilon(self.u)), epsilon(self.v)) * fe.dx
+ fe.dot(self.lamb_, self.u) * fe.dx
+ fe.dot(self.lamb, self.v) * fe.dx
)
self.K = fe.assemble(self.a)
# ! self.solver = fe.LUSolver(self.K) #! Précédemment
# ! Utilisé pour les études de convergence mais pas le plus approprié selon Jérémy
# ! self.solver = fe.KrylovSolver(self.K, method="cg")
self.solver = fe.LUSolver(self.K, "mumps")
self.solver.parameters["symmetric"] = True
# print("Solver parameters : ")
if self._solver["type"] == "Krylov":
self.solver = fe.KrylovSolver(self.K, self._solver["method"])
elif self._solver["type"] == "LU":
self.solver = fe.LUSolver(self.K, self._solver["method"])
self.solver.parameters["symmetric"] = True
if not self._solver["preconditioner"] is None:
self.solver.parameters.preconditioner = self._solver["preconditioner"]
# fe.info(self.solver.parameters, True)
self.localization = {} # dictionary of localization field objects, to be filled up when calling auxiliary problems (lazy evaluation)
self.ConstitutiveTensors = {} # dictionary of homogenized constitutive tensors, to be filled up when calling auxiliary problems and averaging localization fields.
self.localization = dict()
# dictionary of localization field objects,
# will be filled up when calling auxiliary problems (lazy evaluation)
self.ConstitutiveTensors = dict()
# dictionary of homogenized constitutive tensors,
# will be filled up when calling auxiliary problems and
# averaging localization fields.
def homogenizationScheme(self, model):
DictOfLocalizationsU = {}
......@@ -129,189 +138,219 @@ class Fenics2DHomogenization(object):
DictOfLocalizationsEpsilon = {}
DictOfConstitutiveTensors = {}
DictOfLocalizationsU['U'] = self.LocalizationU()['U']
DictOfLocalizationsU["U"] = self.LocalizationU()["U"]
if model == 'E' or model == 'EG' or model == 'EGG':
DictOfLocalizationsU['E'] = self.LocalizationE()['U']
DictOfLocalizationsSigma['E'] = self.LocalizationE()['Sigma']
DictOfLocalizationsEpsilon['E'] = self.LocalizationE()['Epsilon']
if model == "E" or model == "EG" or model == "EGG":
DictOfLocalizationsU["E"] = self.LocalizationE()["U"]
DictOfLocalizationsSigma["E"] = self.LocalizationE()["Sigma"]
DictOfLocalizationsEpsilon["E"] = self.LocalizationE()["Epsilon"]
if model == 'EG':
DictOfLocalizationsU['EGbis'] = self.LocalizationEG()['U']
DictOfLocalizationsSigma['EGbis'] = self.LocalizationEGbis()['Sigma']
DictOfLocalizationsEpsilon['EGbis'] = self.LocalizationEGbis()['Epsilon']
if model == "EG":
DictOfLocalizationsU["EGbis"] = self.LocalizationEG()["U"]
DictOfLocalizationsSigma["EGbis"] = self.LocalizationEGbis()["Sigma"]
DictOfLocalizationsEpsilon["EGbis"] = self.LocalizationEGbis()["Epsilon"]
if model == 'EGG':
DictOfLocalizationsU['EG'] = self.LocalizationEG()['U']
DictOfLocalizationsSigma['EG'] = self.LocalizationEG()['Sigma']
DictOfLocalizationsEpsilon['EG'] = self.LocalizationEG()['Epsilon']
if model == "EGG":
DictOfLocalizationsU["EG"] = self.LocalizationEG()["U"]
DictOfLocalizationsSigma["EG"] = self.LocalizationEG()["Sigma"]
DictOfLocalizationsEpsilon["EG"] = self.LocalizationEG()["Epsilon"]
DictOfLocalizationsU['EGGbis'] = self.LocalizationEGG()['U']
DictOfLocalizationsSigma['EGGbis'] = self.LocalizationEGGbis()['Sigma']
DictOfLocalizationsEpsilon['EGGbis'] = self.LocalizationEGGbis()['Epsilon']
DictOfLocalizationsU["EGGbis"] = self.LocalizationEGG()["U"]
DictOfLocalizationsSigma["EGGbis"] = self.LocalizationEGGbis()["Sigma"]
DictOfLocalizationsEpsilon["EGGbis"] = self.LocalizationEGGbis()["Epsilon"]
Keys = list(DictOfLocalizationsSigma.keys())
for Key1 in Keys:
DictOfConstitutiveTensors[Key1] = {}
for Key2 in Keys:
DictOfConstitutiveTensors[Key1][Key2] = 0.
DictOfConstitutiveTensors[Key1][Key2] = 0.0
nk = len(Keys)
for i in range(nk):
Key1 = Keys[i]
for j in range(i, nk):
Key2 = Keys[j]
C = self.CrossEnergy(Key1, Key2, DictOfLocalizationsSigma[Key1], DictOfLocalizationsEpsilon[Key2])
C = self.CrossEnergy(
Key1,
Key2,
DictOfLocalizationsSigma[Key1],
DictOfLocalizationsEpsilon[Key2],
)
DictOfConstitutiveTensors[Key1][Key2] = C
DictOfConstitutiveTensors[Key2][Key1] = C.T
return DictOfLocalizationsU, DictOfLocalizationsSigma, DictOfLocalizationsEpsilon, DictOfConstitutiveTensors
return (
DictOfLocalizationsU,
DictOfLocalizationsSigma,
DictOfLocalizationsEpsilon,
DictOfConstitutiveTensors,
)
def LocalizationU(self):
try:
out = self.localization['U']
out = self.localization["U"]
except KeyError:
u = [fe.interpolate(fe.Constant((1., 0.)), self.V),
fe.interpolate(fe.Constant((0., 1.)), self.V)]
s = [fe.interpolate(fe.Constant((0., 0., 0.)), self.W),
fe.interpolate(fe.Constant((0., 0., 0.)), self.W)]
out = {'U': u, 'Sigma': s, 'Epsilon': s}
self.localization['U'] = out
u = [
fe.interpolate(fe.Constant((1.0, 0.0)), self.V),
fe.interpolate(fe.Constant((0.0, 1.0)), self.V),
]
s = [
fe.interpolate(fe.Constant((0.0, 0.0, 0.0)), self.W),
fe.interpolate(fe.Constant((0.0, 0.0, 0.0)), self.W),
]
out = {"U": u, "Sigma": s, "Epsilon": s}
self.localization["U"] = out
return out
def LocalizationE(self):
'''
"""
return E localization fields (as function of macro E)
'''
"""
try:
out = self.localization['E']
out = self.localization["E"]
except KeyError:
f = [fe.Constant((0, 0)), fe.Constant((0, 0)), fe.Constant((0, 0))]
Fload = [fe.interpolate(fo, self.V) for fo in f]
epsilon0 = [fe.Constant((1, 0, 0)), fe.Constant((0, 1, 0)), fe.Constant((0, 0, 1))]
epsilon0 = [
fe.Constant((1, 0, 0)),
fe.Constant((0, 1, 0)),
fe.Constant((0, 0, 1)),
]
Epsilon0 = [fe.interpolate(eps, self.W) for eps in epsilon0]
u, sigm, eps = self.genericAuxiliaryProblem(Fload, Epsilon0)
out = {'U': u, 'Sigma': sigm, 'Epsilon': eps}
self.localization['E'] = out
out = {"U": u, "Sigma": sigm, "Epsilon": eps}
self.localization["E"] = out
return out
def LocalizationEbis(self):
'''
"""
return E stress/strain localization fields $u^U\diads\delta$ (as function of macro E)
'''
"""
try:
out = self.localization['Ebis']
out = self.localization["Ebis"]
except KeyError:
epsilon0 = [fe.Constant((1, 0, 0)), fe.Constant((0, 1, 0)), fe.Constant((0, 0, 1))]
epsilon0 = [
fe.Constant((1, 0, 0)),
fe.Constant((0, 1, 0)),
fe.Constant((0, 0, 1)),
]
Epsilon0 = [fe.interpolate(eps, self.W) for eps in epsilon0]
out = {
'Sigma': [sigma(self.rve.C_per, Epsilon0[i]) for i in range(3)],
'Epsilon': Epsilon0}
self.localization['Ebis'] = out
"Sigma": [sigma(self.rve.C_per, Epsilon0[i]) for i in range(3)],
"Epsilon": Epsilon0,
}
self.localization["Ebis"] = out
return out
def LocalizationEG(self):
'''
"""
return GradE localization fields (as function of macro GradE)
'''
"""
try:
out = self.localization['EG']
out = self.localization["EG"]
except KeyError:
self.LocalizationE()
Fload = self.Sigma2Fload(self.localization['E']['Sigma'])
Epsilon0 = self.Displacement2Epsilon0(self.localization['E']['U'])
Fload = self.Sigma2Fload(self.localization["E"]["Sigma"])
Epsilon0 = self.Displacement2Epsilon0(self.localization["E"]["U"])
u, sigm, eps = self.genericAuxiliaryProblem(Fload, Epsilon0)
out = {'U': u, 'Sigma': sigm, 'Epsilon': eps}
self.localization['EG'] = out
out = {"U": u, "Sigma": sigm, "Epsilon": eps}
self.localization["EG"] = out
return out
def LocalizationEGbis(self):
'''
"""
return GradE stress/strain localization fields for the EG model $u^E\diads\delta$ (as function of macro GradE)
'''
"""
try:
out = self.localization['EGbis']
out = self.localization["EGbis"]
except KeyError:
# h = self.lamination.total_thickness
self.LocalizationE()
Epsilon0 = self.Displacement2Epsilon0(self.localization['E']['U'])
Epsilon0 = self.Displacement2Epsilon0(self.localization["E"]["U"])
out = {
'Sigma': [sigma(self.rve.C_per, Epsilon0[i]) for i in range(6)],
'Epsilon': Epsilon0}
self.localization['EGbis'] = out
"Sigma": [sigma(self.rve.C_per, Epsilon0[i]) for i in range(6)],
"Epsilon": Epsilon0,
}
self.localization["EGbis"] = out
return out
def LocalizationEGG(self):
'''
"""
return GradGradE localization fields (as function of macro GradGradE)
>>> a = LaminatedComposite.paganosLaminate([1.0, 1.5, 1.0], [np.pi/3, -np.pi/3, np.pi/3])
>>> Lam = Laminate2D(a)
>>> plt=plotTransverseDistribution(Lam.LocalizationEGG()['U'])
>>> show()
'''
"""
try:
out = self.localization['EGG']
out = self.localization["EGG"]
except KeyError:
self.LocalizationEG()
Fload = self.Sigma2Fload(self.localization['EG']['Sigma'])
Epsilon0 = self.Displacement2Epsilon0(self.localization['EG']['U'])
Fload = self.Sigma2Fload(self.localization["EG"]["Sigma"])
Epsilon0 = self.Displacement2Epsilon0(self.localization["EG"]["U"])
u, sigm, eps = self.genericAuxiliaryProblem(Fload, Epsilon0)
out = {'U': u, 'Sigma': sigm, 'Epsilon': eps}
self.localization['EGG'] = out
out = {"U": u, "Sigma": sigm, "Epsilon": eps}
self.localization["EGG"] = out
return out
def LocalizationEGGbis(self):
'''
"""
return GradGradE stress/strain localization fields for the EGG model $u^EG\diads\delta$ (as function of macro GradGradE)
>>> a = LaminatedComposite.paganosLaminate([1.0, 1.5, 1.0], [np.pi/3, -np.pi/3, np.pi/3])
>>> Lam = Laminate2D(a)
>>> plt=plotTransverseDistribution(Lam.LocalizationEGGbis()['Sigma'])
>>> show()
'''
"""
try:
out = self.localization['EGGbis']
out = self.localization["EGGbis"]
except KeyError:
# h = self.lamination.total_thickness
self.LocalizationEG()
Epsilon0 = self.Displacement2Epsilon0(self.localization['EG']['U'])
Epsilon0 = self.Displacement2Epsilon0(self.localization["EG"]["U"])
out = {
'Sigma': [sigma(self.rve.C_per, Epsilon0[i]) for i in range(12)],
'Epsilon': Epsilon0}
self.localization['EGGbis'] = out
"Sigma": [sigma(self.rve.C_per, Epsilon0[i]) for i in range(12)],
"Epsilon": Epsilon0,
}
self.localization["EGGbis"] = out
return out
def LocalizationSG(self):
'''
"""
return GradSigma localization fields (as function of macro GradSigma)
'''
"""
try:
out = self.localization['SG']
except:
out = self.localization["SG"]
except KeyError:
self.LocalizationE()
C = self.CrossEnergy('E','E',self.localization['E']['Sigma'], self.localization['E']['Epsilon'])
BEps = self.localization['E']['Sigma']
C = self.CrossEnergy(
"E",
"E",
self.localization["E"]["Sigma"],
self.localization["E"]["Epsilon"],
)
BEps = self.localization["E"]["Sigma"]
n = len(BEps)
BSigma = [[BEps[i][j] * C for j in range(n)] for i in range(n)]
# Assembling the constitutive matrix
......@@ -319,11 +358,13 @@ class Fenics2DHomogenization(object):
for i in range(n):
BSigma_j = []
for j in range(n):
B =fe.Function(self.X)
B = fe.Function(self.X)
for k in range(n):
print(C[k, j])
print(BEps[i].sub(k))
B = fe.interpolate(B + BEps[i].sub(k) * fe.Constant(C[k,j]),self.X)
B = fe.interpolate(
B + BEps[i].sub(k) * fe.Constant(C[k, j]), self.X
)
BSigma_j = BSigma_j + [B]
BSigma_i = BSigma_i + [BSigma_j]
......@@ -333,11 +374,10 @@ class Fenics2DHomogenization(object):
u, sigm, eps = self.genericAuxiliaryProblem(Fload, Epsilon0)
out = {'U': u, 'Sigma': sigm, 'Epsilon': eps}
self.localization['EG'] = out
out = {"U": u, "Sigma": sigm, "Epsilon": eps}
self.localization["SG"] = out
return out
def CrossEnergy(self, Key1, Key2, Loc1, Loc2):
"""
Calcul l'energie croisee des tenseurs de localisation Loc1 et Loc2.
......@@ -353,8 +393,10 @@ class Fenics2DHomogenization(object):
C = np.zeros((len_1, len_1))
for k in range(len_1):
for l in range(k, len_1):
C[k, l] = (cross_energy(Loc1[k], Loc2[l], self.rve.mesh)
/ self.rve.rve_area)
C[k, l] = (
cross_energy(Loc1[k], Loc2[l], self.rve.mesh)
/ self.rve.rve_area
)
C[l, k] = C[k, l]
else:
len_1 = len(Loc1)
......@@ -362,8 +404,10 @@ class Fenics2DHomogenization(object):
C = np.zeros((len_1, len_2))
for k in range(len_1):
for l in range(len_2):
C[k, l] = (cross_energy(Loc1[k], Loc2[l], self.rve.mesh)
/ self.rve.rve_area)
C[k, l] = (
cross_energy(Loc1[k], Loc2[l], self.rve.mesh)
/ self.rve.rve_area
)
try:
self.ConstitutiveTensors[Key1][Key2] = C
except KeyError:
......@@ -390,16 +434,17 @@ class Fenics2DHomogenization(object):
E2 = []
for i in range(len(Fload)):
logger.info(f"Progression : load {i+1} / {len(Fload)}")
L = (fe.dot(Fload[i], self.v)
+ fe.inner(-sigma(self.rve.C_per, Epsilon0[i]), epsilon(self.v))
) * fe.dx
L = (
fe.dot(Fload[i], self.v)
+ fe.inner(-sigma(self.rve.C_per, Epsilon0[i]), epsilon(self.v))
) * fe.dx
# u_s = fe.Function(self.V)
res = fe.assemble(L)
# self.bc.apply(res) #TODO à tester. Pas nécessaire pour le moment, la ligne # K,res = fe.assemble_system(self.a,L,self.bc) était commentée.
# self.solver.solve(K, u_s.vector(), res) #* Previous method
self.solver.solve(self.w.vector(), res)
#* More info : https://fenicsproject.org/docs/dolfin/1.5.0/python/programmers-reference/cpp/la/PETScLUSolver.html
# * More info : https://fenicsproject.org/docs/dolfin/1.5.0/python/programmers-reference/cpp/la/PETScLUSolver.html
u_only = fe.interpolate(self.w.sub(0), self.V)
# ? Essayer de faire fe.assign(self.u_only, self.w.sub(0)) ?
......@@ -422,7 +467,11 @@ class Fenics2DHomogenization(object):
f = [fe.Constant((0, 0)), fe.Constant((0, 0)), fe.Constant((0, 0))]
Fload = [fe.interpolate(fo, self.V) for fo in f]
epsilon0 = [fe.Constant((1, 0, 0)), fe.Constant((0, 1, 0)), fe.Constant((0, 0, 1))]
epsilon0 = [
fe.Constant((1, 0, 0)),
fe.Constant((0, 1, 0)),
fe.Constant((0, 0, 1)),
]
Epsilon0 = [fe.interpolate(eps, self.W) for eps in epsilon0]
U = []
......@@ -432,9 +481,9 @@ class Fenics2DHomogenization(object):
for i in range(order):
logger.info(i)
if i > 0:
logger.info('compute load')
Fload = self.Sigma2Fload(Sigma[i-1])
Epsilon0 = self.Displacement2Epsilon0(U[i-1])
logger.info("compute load")
Fload = self.Sigma2Fload(Sigma[i - 1])
Epsilon0 = self.Displacement2Epsilon0(U[i - 1])
u, sigm, eps = self.genericAuxiliaryProblem(Fload, Epsilon0)
U = U + [u]
......@@ -444,36 +493,52 @@ class Fenics2DHomogenization(object):
return U, Sigma, Epsilon
def Sigma2Fload(self, S):
'''genere le chargement volumique associe au localisateur en
contraintes du probleme precedent'''
"""genere le chargement volumique associe au localisateur en
contraintes du probleme precedent"""
Fload = []
S_avout = []
# Oter la valeur moyenne
for i in range(len(S)):
Sj = []
for j in range(int(self.topo_dim*(self.topo_dim+1)/2)):
S_av = fe.assemble(S[i][j]*fe.dx)/self.rve.mat_area
for j in range(int(self.topo_dim * (self.topo_dim + 1) / 2)):
S_av = fe.assemble(S[i][j] * fe.dx) / self.rve.mat_area
S_av = fe.interpolate(fe.Constant(S_av), self.X)
Sj = Sj + [S[i][j] - S_av]
S_avout = S_avout +[Sj]
S_avout = S_avout + [Sj]
# Allocation du chargement
for i in range(len(S)):
Fload = Fload + [fe.as_vector((S_avout[i][0],S_avout[i][2]))]
Fload = Fload + [fe.as_vector((S_avout[i][0], S_avout[i][2]))]
for i in range(len(S)):
Fload = Fload + [fe.as_vector((S_avout[i][2],S_avout[i][1]))]
Fload = Fload + [fe.as_vector((S_avout[i][2], S_avout[i][1]))]
return Fload
def Displacement2Epsilon0(self, U):
'''Converti le localisateur en deplacement en champ de predeformation
a appliquer au probleme auxiliaire suivant'''
"""Converti le localisateur en deplacement en champ de predeformation
a appliquer au probleme auxiliaire suivant"""
Epsilon0 = []
for i in range(len(U)):
Epsilon0 = Epsilon0 + [fe.as_vector((U[i][0], fe.interpolate(fe.Constant(0.),self.X), U[i][1]/fe.sqrt(2)))]
Epsilon0 = Epsilon0 + [