programmer's documentation
Functions
cs_cdovb_diffusion.c File Reference

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"
Include dependency graph for cs_cdovb_diffusion.c:

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_tcs_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_tcs_cdovb_diffusion_builder_free (cs_cdovb_diff_t *diff)
 Free a cs_cdovb_diff_t structure. More...
 
cs_locmat_tcs_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_tcs_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...
 

Detailed Description

Build discrete stiffness matrices and handled boundary conditions diffusion term in CDO vertex-based schemes.

Function Documentation

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 
)
static

Compute the stiffness matrix on this cell using a Whitney Barycentric Subdivision (WBS) algo.

Parameters
[in]c_idcurrent cell id
[in]connectpointer to a cs_cdo_connect_t struct.
[in]quantpointer to a cs_cdo_quantities_t struct.
[in]loc_idslocal vertex numbering on a cell
[in]tensor3x3 matrix attached to the diffusion property
[in,out]diffauxiliary structure used to build the diff. term
Returns
a pointer to a local stiffness matrix
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 
)
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.

Parameters
[in]f_idid of the face
[in]pfqgeometric quantities related to face f_id
[in]deqgeometric quantities related to the dual edge
[in]connectpointer to a cs_cdo_connect_t structure
[in]quantpointer to a cs_cdo_quantites_t structure
[in]loc_idslocal vertex numbering on a cell
[in,out]wvfpointer to an array storing the weight/vertex
[in,out]pefc_volpointer to an array storing the volume of pefc
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 
)
static

Compute and store the unit vector tangential to x_c –> x_v. Store also the local weight related to each vertex wvc.

Parameters
[in]c_idcurrent cell id
[in]connectpointer to a cs_cdo_connect_t struct.
[in]quantpointer to a cs_cdo_quantities_t struct.
[in,out]diffauxiliary structure used to build the diff. term
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 
)
static

Encode E_c cap E_v sets for all vertices of a cell using a mask of bits.

(end ignore by Doxygen)

Parameters
[in]c_idcell id
[in]connectpoint to a cs_cdo_connect_t struct.
[in]vtagtag on vertices to store local ids
[in,out]diffpointer to a cs_cdovb_diff_t structure
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 
)
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.

Parameters
[in]ifcincidence number i_{f,c}
[in]pfqprimal face quantities
[in]deqdual edge quantities
[in,out]grd_funcgradient of the Lagrange function related to c
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 
)
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'.

Parameters
[in]ubaseunit vector from xc to xv' (opposite)
[in]udirunit vector from xc to xv
[in]len_vclenght of the segment [xc, xv]
[in]deqdual edge quantities
[in,out]grd_funcgradient 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.

Parameters
[in]c_idcurrent cell id
[in]connectpointer to a cs_cdo_connect_t struct.
[in]quantpointer to a cs_cdo_quantities_t struct.
[in]vtagpointer to a cs_cdovb_scaleq_t struct.
[in]tensor3x3 matrix attached to the diffusion property
[in,out]diffauxiliary structure used to build the diff. term
Returns
a pointer to a local stiffness matrix
cs_cdovb_diff_t* cs_cdovb_diffusion_builder_free ( cs_cdovb_diff_t diff)

Free a cs_cdovb_diff_t structure.

Parameters
[in,out]diffpointer to a cs_cdovb_diff_t struc.
Returns
NULL
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.

Parameters
[in]connectpointer to a cs_cdo_connect_t struct.
[in]is_uniformdiffusion tensor is uniform ? (true or false)
[in]h_infocs_param_hodge_t struct.
[in]bc_enforcetype of boundary enforcement for Dirichlet values
Returns
a pointer to a new allocated cs_cdovb_diffusion_builder_t struc.
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.

Parameters
[in]c_idcell id
[in]f_idface id (a border face attached to a Dir. BC)
[in]connectpointer to a cs_cdo_connect_t struct.
[in]quantpointer to a cs_cdo_quantities_t struct.
[in]matpty3x3 matrix related to the diffusion property
[in]eig_ratioeigenvalue_max/eigenvalue_min
[in]eig_maxeigenvalue with maximal value
[in,out]loc_v_idsstore local vertex ids
[in,out]v_coefstore local contribution on each border vertex
[in,out]diffauxiliary structure used to build the diff. term
Returns
a pointer to a local "normal trace gradient" matrix