jax_fem.fe
module#
- class jax_fem.fe.FiniteElement(mesh, vec, dim, ele_type, gauss_order, dirichlet_bc_info)[source]#
Finite element class for one variable. This variable can be:
- Parameters:
mesh (Mesh)
vec (int)
dim (int)
ele_type (str)
gauss_order (int)
dirichlet_bc_info (list)
- vec#
Number of vector components in solution (e.g., 3 for displacement in 3D, 1 for temperature).
- Type:
int
- dim#
Spatial dimension of the problem. Currently supports:
2 (2D problems)
3 (3D problems)
- Type:
int
- ele_type#
Element type. Currently supports:
‘QUAD4’
‘QUAD8’
‘TRI3’
‘TRI6’
‘HEX8’
‘HEX20’
‘HEX27’
‘TET4’
‘TET10’
- Type:
str
- gauss_order#
Order of Gaussian quadrature.
- Type:
int
- dirichlet_bc_info#
A list for Dirichlet boundary condition information, whose elements are structured as:
location_fns: list of callables Each callable takes a point (NumpyArray) and returns a boolean indicating if the point satisfies the location condition
vecs: list of integers Each integer must be in the range of 0 to vec - 1, specifying which component of the (vector) variable to apply Dirichlet condition to
value_fns: list of callables Each callable takes a point and returns the Dirichlet value
- Type:
list
- get_shape_grads()[source]#
Compute shape function gradient value.
The gradient is w.r.t physical coordinates. Refer to Hughes, Thomas JR. The finite element method: linear static and dynamic finite element analysis. Courier Corporation, 2012. Page 147, Eq. (3.9.3)
- Returns:
shape_grads_physical (NumpyArray) – Shape is (num_cells, num_quads, num_nodes, dim).
JxW (NumpyArray) – Shape is (num_cells, num_quads).
- get_face_shape_grads(boundary_inds)[source]#
Face shape function gradients and JxW (for surface integral). Nanson’s formula is used to map physical surface ingetral to reference domain. Refer to wikiversity.
- Parameters:
boundary_inds (list[NumpyArray]) – Shape is (num_selected_faces, 2).
- Returns:
face_shape_grads_physical (NumpyArray) – Shape is (num_selected_faces, num_face_quads, num_nodes, dim).
nanson_scale (NumpyArray) – Shape is (num_selected_faces, num_face_quads).
- get_physical_quad_points()[source]#
Compute physical quadrature points
- Returns:
physical_quad_points – Shape is (num_cells, num_quads, dim).
- Return type:
NumpyArray
- get_physical_surface_quad_points(boundary_inds)[source]#
Compute physical quadrature points on the surface
- Parameters:
boundary_inds (list[NumpyArray]) – Shape for NumpyArray is (num_selected_faces, 2).
- Returns:
physical_surface_quad_points – Shape for NumpyArray is (num_selected_faces, num_face_quads, dim).
- Return type:
NumpyArray
- Dirichlet_boundary_conditions(dirichlet_bc_info)[source]#
Indices and values for Dirichlet B.C.
- Parameters:
dirichlet_bc_info (list) – [location_fns, vecs, value_fns]
- Returns:
node_inds_list (list[NumpyArray]) – The value of the NumpyArray ranges from 0 to num_total_nodes - 1.
vec_inds_list (list[NumpyArray]) – The value of the NumpyArray ranges from 0 to to vec - 1.
vals_list (list[NumpyArray]) – Dirichlet values to be assigned.
- update_Dirichlet_boundary_conditions(dirichlet_bc_info)[source]#
Reset Dirichlet boundary conditions. Useful when a time-dependent problem is solved, and at each iteration the boundary condition needs to be updated.
- Parameters:
dirichlet_bc_info (list) – [location_fns, vecs, value_fns]
- get_boundary_conditions_inds(location_fns)[source]#
Given location functions, compute which faces satisfy the condition.
- Parameters:
location_fns (list[callable]) –
The callable is a location function that inputs a point (NumpyArray) and returns if the point satisfies the location condition. For example,
lambda x: np.isclose(x[0], 0.)
If this location function takes 2 arguments, then the first is point and the second is index. For example,
lambda x, ind: np.isclose(x[0], 0.) & np.isin(ind, np.array([1, 3, 10]))
- Returns:
boundary_inds_list – Shape of NumpyArray is (num_selected_faces, 2).
boundary_inds_list[k][i, 0] returns the global cell index of the ith selected face of boundary subset k.
boundary_inds_list[k][i, 1] returns the local face index of the ith selected face of boundary subset k.
- Return type:
list[NumpyArray]
- convert_from_dof_to_quad(sol)[source]#
Obtain quad values from nodal solution.
- Parameters:
sol (JaxArray) – Shape is (num_total_nodes, vec).
- Returns:
u – Shape is (num_cells, num_quads, vec).
- Return type:
JaxArray