![]() |
programmer's documentation
|
Build discrete stiffness matrices and handled boundary conditions diffusion term in CDO vertex-based schemes. More...
#include "cs_defs.h"
#include <stdio.h>
#include <stdlib.h>
#include <assert.h>
#include <math.h>
#include <limits.h>
#include <bft_mem.h>
#include <bft_printf.h>
#include "cs_hodge.h"
#include "cs_property.h"
#include "cs_cdovb_diffusion.h"
Functions | |
static void | _encode_edge_masks (cs_lnum_t c_id, const cs_cdo_connect_t *connect, const cs_lnum_t vtag[], cs_cdovb_diff_t *diff) |
Encode E_c cap E_v sets for all vertices of a cell using a mask of bits. More... | |
static void | _define_wbs_builder (cs_lnum_t c_id, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, cs_cdovb_diff_t *diff) |
Compute and store the unit vector tangential to x_c –> x_v. Store also the local weight related to each vertex wvc. More... | |
static void | _grad_lagrange_cell_pfc (short int ifc, const cs_quant_t pfq, const cs_nvec3_t deq, cs_real_3_t grd_func) |
Compute the gradient of a Lagrange hat function related to a cell c in a p_{f, c} subvolume of a cell where f is a face of c. More... | |
static void | _grad_lagrange_vtx_pefc (const cs_real_3_t ubase, const cs_real_3_t udir, cs_real_t len_vc, const cs_nvec3_t deq, cs_real_3_t grd_func) |
Compute the gradient of a Lagrange hat function related to primal vertices in a p_{e,f, c} subvolume of a cell c where e is an edge belonging to the face f with vertices v and v'. More... | |
static void | _compute_wvf_pefcvol (cs_lnum_t f_id, const cs_quant_t pfq, const cs_nvec3_t deq, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_lnum_t *loc_ids, cs_real_t *wvf, cs_real_t *pefc_vol) |
Compute for each face a weight related to each vertex w_{v,f} This weight is equal to |dc(v) cap f|/|f| so that the sum of the weights is equal to 1. Set also the local and local numbering of the vertices of this face. More... | |
static cs_locmat_t * | _compute_wbs_stiffness (cs_lnum_t c_id, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_lnum_t *loc_ids, const cs_real_3_t *tensor, cs_cdovb_diff_t *diff) |
Compute the stiffness matrix on this cell using a Whitney Barycentric Subdivision (WBS) algo. More... | |
cs_cdovb_diff_t * | cs_cdovb_diffusion_builder_init (const cs_cdo_connect_t *connect, bool is_uniform, const cs_param_hodge_t h_info, const cs_param_bc_enforce_t bc_enforce) |
Initialize a builder structure used to build the stiffness matrix. More... | |
cs_cdovb_diff_t * | cs_cdovb_diffusion_builder_free (cs_cdovb_diff_t *diff) |
Free a cs_cdovb_diff_t structure. More... | |
cs_locmat_t * | cs_cdovb_diffusion_build_local (cs_lnum_t c_id, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_lnum_t *vtag, const cs_real_3_t *tensor, cs_cdovb_diff_t *diff) |
Define the local (cellwise) stiffness matrix. More... | |
cs_locmat_t * | cs_cdovb_diffusion_ntrgrd_build (cs_lnum_t c_id, cs_lnum_t f_id, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_real_t matpty[3][3], cs_real_t eig_ratio, cs_real_t eig_max, cs_lnum_t *loc_v_ids, cs_real_t *v_coef, cs_cdovb_diff_t *diff) |
Define the local (cellwise) "normal trace gradient" matrix This local matrix is used in Nitsche method to weakly penalized Dirichlet boundary conditions. More... | |
Build discrete stiffness matrices and handled boundary conditions diffusion term in CDO vertex-based schemes.
|
static |
Compute the stiffness matrix on this cell using a Whitney Barycentric Subdivision (WBS) algo.
[in] | c_id | current cell id |
[in] | connect | pointer to a cs_cdo_connect_t struct. |
[in] | quant | pointer to a cs_cdo_quantities_t struct. |
[in] | loc_ids | local vertex numbering on a cell |
[in] | tensor | 3x3 matrix attached to the diffusion property |
[in,out] | diff | auxiliary structure used to build the diff. term |
|
static |
Compute for each face a weight related to each vertex w_{v,f} This weight is equal to |dc(v) cap f|/|f| so that the sum of the weights is equal to 1. Set also the local and local numbering of the vertices of this face.
[in] | f_id | id of the face |
[in] | pfq | geometric quantities related to face f_id |
[in] | deq | geometric quantities related to the dual edge |
[in] | connect | pointer to a cs_cdo_connect_t structure |
[in] | quant | pointer to a cs_cdo_quantites_t structure |
[in] | loc_ids | local vertex numbering on a cell |
[in,out] | wvf | pointer to an array storing the weight/vertex |
[in,out] | pefc_vol | pointer to an array storing the volume of pefc |
|
static |
Compute and store the unit vector tangential to x_c –> x_v. Store also the local weight related to each vertex wvc.
[in] | c_id | current cell id |
[in] | connect | pointer to a cs_cdo_connect_t struct. |
[in] | quant | pointer to a cs_cdo_quantities_t struct. |
[in,out] | diff | auxiliary structure used to build the diff. term |
|
static |
Encode E_c cap E_v sets for all vertices of a cell using a mask of bits.
(end ignore by Doxygen)
[in] | c_id | cell id |
[in] | connect | point to a cs_cdo_connect_t struct. |
[in] | vtag | tag on vertices to store local ids |
[in,out] | diff | pointer to a cs_cdovb_diff_t structure |
|
static |
Compute the gradient of a Lagrange hat function related to a cell c in a p_{f, c} subvolume of a cell where f is a face of c.
[in] | ifc | incidence number i_{f,c} |
[in] | pfq | primal face quantities |
[in] | deq | dual edge quantities |
[in,out] | grd_func | gradient of the Lagrange function related to c |
|
static |
Compute the gradient of a Lagrange hat function related to primal vertices in a p_{e,f, c} subvolume of a cell c where e is an edge belonging to the face f with vertices v and v'.
[in] | ubase | unit vector from xc to xv' (opposite) |
[in] | udir | unit vector from xc to xv |
[in] | len_vc | lenght of the segment [xc, xv] |
[in] | deq | dual edge quantities |
[in,out] | grd_func | gradient of the Lagrange shape function |
cs_locmat_t* cs_cdovb_diffusion_build_local | ( | cs_lnum_t | c_id, |
const cs_cdo_connect_t * | connect, | ||
const cs_cdo_quantities_t * | quant, | ||
const cs_lnum_t * | vtag, | ||
const cs_real_3_t * | tensor, | ||
cs_cdovb_diff_t * | diff | ||
) |
Define the local (cellwise) stiffness matrix.
[in] | c_id | current cell id |
[in] | connect | pointer to a cs_cdo_connect_t struct. |
[in] | quant | pointer to a cs_cdo_quantities_t struct. |
[in] | vtag | pointer to a cs_cdovb_scaleq_t struct. |
[in] | tensor | 3x3 matrix attached to the diffusion property |
[in,out] | diff | auxiliary structure used to build the diff. term |
cs_cdovb_diff_t* cs_cdovb_diffusion_builder_free | ( | cs_cdovb_diff_t * | diff | ) |
Free a cs_cdovb_diff_t structure.
[in,out] | diff | pointer to a cs_cdovb_diff_t struc. |
cs_cdovb_diff_t* cs_cdovb_diffusion_builder_init | ( | const cs_cdo_connect_t * | connect, |
bool | is_uniform, | ||
const cs_param_hodge_t | h_info, | ||
const cs_param_bc_enforce_t | bc_enforce | ||
) |
Initialize a builder structure used to build the stiffness matrix.
[in] | connect | pointer to a cs_cdo_connect_t struct. |
[in] | is_uniform | diffusion tensor is uniform ? (true or false) |
[in] | h_info | cs_param_hodge_t struct. |
[in] | bc_enforce | type of boundary enforcement for Dirichlet values |
cs_locmat_t* cs_cdovb_diffusion_ntrgrd_build | ( | cs_lnum_t | c_id, |
cs_lnum_t | f_id, | ||
const cs_cdo_connect_t * | connect, | ||
const cs_cdo_quantities_t * | quant, | ||
const cs_real_t | matpty[3][3], | ||
cs_real_t | eig_ratio, | ||
cs_real_t | eig_max, | ||
cs_lnum_t * | loc_v_ids, | ||
cs_real_t * | v_coef, | ||
cs_cdovb_diff_t * | diff | ||
) |
Define the local (cellwise) "normal trace gradient" matrix This local matrix is used in Nitsche method to weakly penalized Dirichlet boundary conditions.
[in] | c_id | cell id |
[in] | f_id | face id (a border face attached to a Dir. BC) |
[in] | connect | pointer to a cs_cdo_connect_t struct. |
[in] | quant | pointer to a cs_cdo_quantities_t struct. |
[in] | matpty | 3x3 matrix related to the diffusion property |
[in] | eig_ratio | eigenvalue_max/eigenvalue_min |
[in] | eig_max | eigenvalue with maximal value |
[in,out] | loc_v_ids | store local vertex ids |
[in,out] | v_coef | store local contribution on each border vertex |
[in,out] | diff | auxiliary structure used to build the diff. term |