-
Notifications
You must be signed in to change notification settings - Fork 1
/
my_classes.py
178 lines (137 loc) · 5.42 KB
/
my_classes.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
from mpi4py import MPI
from petsc4py import PETSc
from dolfinx.io import XDMFFile
from dolfinx.fem import (
Constant,
)
from ufl import (
Measure,
)
class Mesh:
"""
Mesh class
Attributes:
mesh (dolfinx.Mesh): the mesh
volume_markers (dolfinx.MeshTags): markers of the mesh cells
surface_markers (dolfinx.MeshTags): markers of the mesh facets
dx (dolfinx.Measure):
ds (dolfinx.Measure):
"""
def __init__(
self, mesh=None, volume_markers=None, surface_markers=None, subdomains=[]
) -> None:
"""Inits Mesh
Args:
mesh (dolfinx.Mesh, optional): the mesh. Defaults to None.
volume_markers (dolfinx.MeshTags, optional): markers of the mesh cells. Defaults to None.
surface_markers (dolfinx.MeshTags, optional): markers of the mesh facets. Defaults to None.
subdomains (list, optional): list of festimx.Subdomain objects
"""
self.mesh = mesh
self.volume_markers = volume_markers
self.surface_markers = surface_markers
self.subdomains = subdomains
self.dx = None
self.ds = None
if self.mesh is not None:
# create cell to facet connectivity
self.mesh.topology.create_connectivity(
self.mesh.topology.dim, self.mesh.topology.dim - 1
)
# create facet to cell connectivity
self.mesh.topology.create_connectivity(
self.mesh.topology.dim - 1, self.mesh.topology.dim
)
def define_markers(self):
self.volume_markers = self.define_volume_markers()
self.surface_markers = self.define_surface_markers()
def define_measures(self):
"""Creates the ufl.Measure objects for self.dx and self.ds"""
self.ds = Measure("ds", domain=self.mesh, subdomain_data=self.surface_markers)
self.dx = Measure("dx", domain=self.mesh, subdomain_data=self.volume_markers)
class MeshXDMF(Mesh):
def __init__(self, cell_file, facet_file, subdomains) -> None:
"""_summary_
Args:
cell_file (str): _description_
facet_file (str): _description_
subdomains (list, optional): list of festimx.Subdomain objects
"""
self.cell_file = cell_file
self.facet_file = facet_file
volume_file = XDMFFile(MPI.COMM_WORLD, self.cell_file, "r")
mesh = volume_file.read_mesh(name="Grid")
super().__init__(mesh=mesh, subdomains=subdomains)
def define_surface_markers(self):
"""Creates the surface markers
Returns:
dolfinx.MeshTags: the tags containing the surface
markers
"""
bondary_file = XDMFFile(MPI.COMM_WORLD, self.facet_file, "r")
mesh_tags_facets = bondary_file.read_meshtags(self.mesh, name="Grid")
return mesh_tags_facets
def define_volume_markers(self):
"""Creates the volume markers
Returns:
dolfinx.MeshTags: the tags containing the volume
markers
"""
volume_file = XDMFFile(MPI.COMM_WORLD, self.cell_file, "r")
mesh_tags_cells = volume_file.read_meshtags(self.mesh, name="Grid")
return mesh_tags_cells
class BoundaryCondition:
"""Base boundary condition class
Args:
surfaces (list, int): the surfaces of the BC
"""
def __init__(self, surfaces) -> None:
if not isinstance(surfaces, list):
surfaces = [surfaces]
self.surfaces = surfaces
class FluxBC(BoundaryCondition):
"""Boundary condition applying flux
Args:
surfaces (list, int): surfaces of the BC
value (float, int, ufl.expression): value of the flux. Deaults to
None
"""
def __init__(self, surfaces, value=None, **kwargs) -> None:
super().__init__(surfaces=surfaces, **kwargs)
self.value = value
def create_form(self, T):
self.form = self.value
class ConvenctiveFluxBC(FluxBC):
"""fluxbc subclass for convective heat flux
-lambda * grad(T) * n = h_coeff * (T - T_ext)
Args:
h_coeff (float, ufl.expression): heat transfer coefficient (W/ms/K)
T_ext (float, ufl.expression): fluid bulk temperature (K)
surfaces (list, int): surfaces of the BC
"""
def __init__(self, h_coeff, T_ext, surfaces) -> None:
super().__init__(surfaces=surfaces)
self.h_coeff = h_coeff
self.T_ext = T_ext
def create_form(self, T):
self.form = -self.h_coeff * (T - self.T_ext)
class Source:
"""Volumetric source term
Args:
value: (float, int, ufl.expression): the value of the
volumteric source term
volume (int, list): the volume in which the source is
applied
mesh (fenics.mesh, optional): the domain in which the
constant value is applied
"""
def __init__(self, value, volumes, mesh=None) -> None:
if not isinstance(volumes, list):
volumes = [volumes]
self.volumes = volumes
if isinstance(value, (float, int)):
if mesh is None:
raise ValueError("Mesh domain needs to be defined")
self.value = Constant(mesh, PETSc.ScalarType(value))
else:
self.value = value