8.0
general documentation
Loading...
Searching...
No Matches
cs_matrix_building.c File Reference
#include "cs_defs.h"
#include <assert.h>
#include <errno.h>
#include <stdio.h>
#include <stdarg.h>
#include <string.h>
#include <math.h>
#include <float.h>
#include "bft_error.h"
#include "bft_mem.h"
#include "bft_printf.h"
#include "cs_blas.h"
#include "cs_halo.h"
#include "cs_halo_perio.h"
#include "cs_log.h"
#include "cs_mesh.h"
#include "cs_field.h"
#include "cs_gradient.h"
#include "cs_ext_neighborhood.h"
#include "cs_mesh_quantities.h"
#include "cs_parameters.h"
#include "cs_porous_model.h"
#include "cs_prototypes.h"
#include "cs_timer.h"
#include "cs_parall.h"
#include "cs_matrix_building.h"
Include dependency graph for cs_matrix_building.c:

Functions

void matrix (const int *iconvp, const int *idiffp, const int *ndircp, const int *isym, const cs_real_t *thetap, const int *imucpp, const cs_real_t coefbp[], const cs_real_t cofbfp[], const cs_real_t rovsdt[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t xcpp[], cs_real_t da[], cs_real_t xa[])
void matrdt (const int *const iconvp, const int *const idiffp, const int *const isym, const cs_real_t coefbp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t da[])
void cs_matrix_wrapper_scalar (int iconvp, int idiffp, int ndircp, int isym, double thetap, int imucpp, const cs_real_t coefbp[], const cs_real_t cofbfp[], const cs_real_t rovsdt[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t xcpp[], cs_real_t da[], cs_real_t xa[])
void cs_matrix_wrapper_vector (int iconvp, int idiffp, int tensorial_diffusion, int ndircp, int isym, cs_lnum_t eb_size, double thetap, const cs_real_33_t coefbp[], const cs_real_33_t cofbfp[], const cs_real_33_t fimp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_33_t da[], cs_real_t xa[])
void cs_matrix_wrapper_tensor (int iconvp, int idiffp, int tensorial_diffusion, int ndircp, int isym, double thetap, const cs_real_66_t coefbts[], const cs_real_66_t cofbfts[], const cs_real_66_t fimp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_66_t da[], cs_real_t xa[])
void cs_sym_matrix_scalar (const cs_mesh_t *m, int idiffp, double thetap, const cs_real_t cofbfp[], const cs_real_t rovsdt[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t *restrict da, cs_real_t *restrict xa)
 Build the diffusion matrix for a scalar field. (symmetric matrix).
void cs_matrix_scalar (const cs_mesh_t *m, int iconvp, int idiffp, double thetap, int imucpp, const cs_real_t coefbp[], const cs_real_t cofbfp[], const cs_real_t rovsdt[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t xcpp[], cs_real_t *restrict da, cs_real_2_t *restrict xa)
 Build the advection/diffusion matrix for a scalar field (non-symmetric matrix).
void cs_sym_matrix_vector (const cs_mesh_t *m, int idiffp, double thetap, const cs_real_33_t cofbfp[], const cs_real_33_t fimp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_33_t *restrict da, cs_real_t *restrict xa)
 Build the diffusion matrix for a vector field (symmetric matrix).
void cs_sym_matrix_tensor (const cs_mesh_t *m, int idiffp, double thetap, const cs_real_66_t cofbfts[], const cs_real_66_t fimp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_66_t *restrict da, cs_real_t *restrict xa)
 Build the diffusion matrix for a tensor field (symmetric matrix).
void cs_matrix_vector (const cs_mesh_t *m, const cs_mesh_quantities_t *mq, int iconvp, int idiffp, cs_lnum_t eb_size, double thetap, const cs_real_33_t coefbp[], const cs_real_33_t cofbfp[], const cs_real_33_t fimp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_33_t *restrict da, cs_real_2_t *restrict xa)
 Build the advection/diffusion matrix for a vector field (non-symmetric matrix).
void cs_matrix_tensor (const cs_mesh_t *m, int iconvp, int idiffp, double thetap, const cs_real_66_t coefbts[], const cs_real_66_t cofbfts[], const cs_real_66_t fimp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_66_t *restrict da, cs_real_2_t *restrict xa)
 Build the advection/diffusion matrix for a tensor field (non-symmetric matrix).
void cs_matrix_time_step (const cs_mesh_t *m, int iconvp, int idiffp, int isym, const cs_real_t coefbp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t *restrict da)
 Build the diagonal of the advection/diffusion matrix for determining the variable time step, flow, Fourier.
void cs_matrix_anisotropic_diffusion (const cs_mesh_t *m, const cs_mesh_quantities_t *mq, int iconvp, int idiffp, double thetap, const cs_real_33_t coefbp[], const cs_real_33_t cofbfp[], const cs_real_33_t fimp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_33_t i_visc[], const cs_real_t b_visc[], cs_real_33_t *restrict da, cs_real_332_t *restrict xa)
 Build the advection/diffusion matrix for a vector field with a tensorial diffusivity.
void cs_matrix_anisotropic_diffusion_tensor (const cs_mesh_t *m, int iconvp, int idiffp, double thetap, const cs_real_66_t coefbp[], const cs_real_66_t cofbfp[], const cs_real_66_t fimp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_66_t i_visc[], const cs_real_t b_visc[], cs_real_66_t *restrict da, cs_real_662_t *restrict xa)
 Build the advection/diffusion matrix for a tensor field with a tensorial diffusivity.
void cs_sym_matrix_anisotropic_diffusion (const cs_mesh_t *m, int idiffp, double thetap, const cs_real_33_t cofbfp[], const cs_real_33_t fimp[], const cs_real_33_t i_visc[], const cs_real_t b_visc[], cs_real_33_t *restrict da, cs_real_33_t *restrict xa)
 Build the diffusion matrix for a vector field with a tensorial diffusivity (symmetric matrix).
void cs_sym_matrix_anisotropic_diffusion_tensor (const cs_mesh_t *m, int idiffp, double thetap, const cs_real_66_t cofbfp[], const cs_real_66_t fimp[], const cs_real_66_t i_visc[], const cs_real_t b_visc[], cs_real_66_t *restrict da, cs_real_66_t *restrict xa)
 Build the diffusion matrix for a tensor field with a tensorial diffusivity (symmetric matrix).

Function Documentation

◆ cs_matrix_anisotropic_diffusion()

void cs_matrix_anisotropic_diffusion ( const cs_mesh_t * m,
const cs_mesh_quantities_t * mq,
int iconvp,
int idiffp,
double thetap,
const cs_real_33_t coefbp[],
const cs_real_33_t cofbfp[],
const cs_real_33_t fimp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_33_t i_visc[],
const cs_real_t b_visc[],
cs_real_33_t *restrict da,
cs_real_332_t *restrict xa )

Build the advection/diffusion matrix for a vector field with a tensorial diffusivity.

The advection is upwind, the diffusion is not reconstructed. The matrix is split into a diagonal block (3x3 times number of cells) and an extra diagonal part (of dimension 2 times 3x3 the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]mqpointer to mesh quantities structure
[in]iconvpindicator
  • 1 advection
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]coefbpboundary condition array for the variable (implicit part - 3x3 tensor array)
[in]cofbfpboundary condition array for the variable flux (implicit part - 3x3 tensor array)
[in]fimppart of the diagonal
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at border faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_matrix_anisotropic_diffusion_tensor()

void cs_matrix_anisotropic_diffusion_tensor ( const cs_mesh_t * m,
int iconvp,
int idiffp,
double thetap,
const cs_real_66_t coefbp[],
const cs_real_66_t cofbfp[],
const cs_real_66_t fimp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_66_t i_visc[],
const cs_real_t b_visc[],
cs_real_66_t *restrict da,
cs_real_662_t *restrict xa )

Build the advection/diffusion matrix for a tensor field with a tensorial diffusivity.

The advection is upwind, the diffusion is not reconstructed. The matrix is split into a diagonal block (3x3 times number of cells) and an extra diagonal part (of dimension 2 times 3x3 the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]iconvpindicator
  • 1 advection
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]coefbpboundary condition array for the variable (implicit part - 3x3 tensor array)
[in]cofbfpboundary condition array for the variable flux (implicit part - 3x3 tensor array)
[in]fimppart of the diagonal
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at border faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_matrix_scalar()

void cs_matrix_scalar ( const cs_mesh_t * m,
int iconvp,
int idiffp,
double thetap,
int imucpp,
const cs_real_t coefbp[],
const cs_real_t cofbfp[],
const cs_real_t rovsdt[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
const cs_real_t xcpp[],
cs_real_t *restrict da,
cs_real_2_t *restrict xa )

Build the advection/diffusion matrix for a scalar field (non-symmetric matrix).

Build the advection/diffusion matrix for a scalar field.

The advection is upwind, the diffusion is not reconstructed. The matrix is split into a diagonal block (number of cells) and an extra diagonal part (of dimension 2 time the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]iconvpindicator
  • 1 advection
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]imucppindicator
  • 0 do not multiply the convective term by Cp
  • 1 do multiply the convective term by Cp
[in]coefbpboundary condition array for the variable (implicit part)
[in]cofbfpboundary condition array for the variable flux (implicit part)
[in]rovsdtworking array
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at border faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ S_\fib $ at border faces for the matrix
[in]xcpparray of specific heat (Cp)
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_matrix_tensor()

void cs_matrix_tensor ( const cs_mesh_t * m,
int iconvp,
int idiffp,
double thetap,
const cs_real_66_t coefbts[],
const cs_real_66_t cofbfts[],
const cs_real_66_t fimp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_66_t *restrict da,
cs_real_2_t *restrict xa )

Build the advection/diffusion matrix for a tensor field (non-symmetric matrix).

The advection is upwind, the diffusion is not reconstructed. The matrix is split into a diagonal block (6x6 times number of cells) and an extra diagonal part (of dimension 2 time the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]iconvpindicator
  • 1 advection
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]coefbtsboundary condition array for the variable (Implicit part - 6x6 tensor array)
[in]cofbftsboundary condition array for the variable flux (Implicit part - 6x6 tensor array)
[in]fimppart of the diagonal
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at border faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_matrix_time_step()

void cs_matrix_time_step ( const cs_mesh_t * m,
int iconvp,
int idiffp,
int isym,
const cs_real_t coefbp[],
const cs_real_t cofbfp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_t *restrict da )

Build the diagonal of the advection/diffusion matrix for determining the variable time step, flow, Fourier.

Parameters
[in]mpointer to mesh structure
[in]iconvpindicator
  • 1 advection
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]isymindicator
  • 1 symmetric matrix
  • 2 non symmmetric matrix
[in]coefbpboundary condition array for the variable (implicit part)
[in]cofbfpboundary condition array for the variable flux (implicit part)
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at border faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ S_\fib $ at border faces for the matrix
[out]dadiagonal part of the matrix

◆ cs_matrix_vector()

void cs_matrix_vector ( const cs_mesh_t * m,
const cs_mesh_quantities_t * mq,
int iconvp,
int idiffp,
cs_lnum_t eb_size,
double thetap,
const cs_real_33_t coefbp[],
const cs_real_33_t cofbfp[],
const cs_real_33_t fimp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_33_t *restrict da,
cs_real_2_t *restrict xa )

Build the advection/diffusion matrix for a vector field (non-symmetric matrix).

The advection is upwind, the diffusion is not reconstructed. The matrix is split into a diagonal block (3x3 times number of cells) and an extra diagonal part (of dimension 2 time the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]mqpointer to mesh quantities structure
[in]iconvpindicator
  • 1 advection
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]coefbpboundary condition array for the variable (implicit part - 3x3 tensor array)
[in]cofbfpboundary condition array for the variable flux (implicit part - 3x3 tensor array)
[in]fimppart of the diagonal
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at border faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_matrix_wrapper_scalar()

void cs_matrix_wrapper_scalar ( int iconvp,
int idiffp,
int ndircp,
int isym,
double thetap,
int imucpp,
const cs_real_t coefbp[],
const cs_real_t cofbfp[],
const cs_real_t rovsdt[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
const cs_real_t xcpp[],
cs_real_t da[],
cs_real_t xa[] )

◆ cs_matrix_wrapper_tensor()

void cs_matrix_wrapper_tensor ( int iconvp,
int idiffp,
int tensorial_diffusion,
int ndircp,
int isym,
double thetap,
const cs_real_66_t coefbts[],
const cs_real_66_t cofbfts[],
const cs_real_66_t fimp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_66_t da[],
cs_real_t xa[] )

◆ cs_matrix_wrapper_vector()

void cs_matrix_wrapper_vector ( int iconvp,
int idiffp,
int tensorial_diffusion,
int ndircp,
int isym,
cs_lnum_t eb_size,
double thetap,
const cs_real_33_t coefbp[],
const cs_real_33_t cofbfp[],
const cs_real_33_t fimp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_33_t da[],
cs_real_t xa[] )

◆ cs_sym_matrix_anisotropic_diffusion()

void cs_sym_matrix_anisotropic_diffusion ( const cs_mesh_t * m,
int idiffp,
double thetap,
const cs_real_33_t cofbfp[],
const cs_real_33_t fimp[],
const cs_real_33_t i_visc[],
const cs_real_t b_visc[],
cs_real_33_t *restrict da,
cs_real_33_t *restrict xa )

Build the diffusion matrix for a vector field with a tensorial diffusivity (symmetric matrix).

The diffusion is not reconstructed. The matrix is split into a diagonal block (3x3 times number of cells) and an extra diagonal part (of dimension 3x3 the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]cofbfpboundary condition array for the variable flux (implicit part - 3x3 tensor array)
[in]fimp$ \tens{f_s}^{imp} $
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_sym_matrix_anisotropic_diffusion_tensor()

void cs_sym_matrix_anisotropic_diffusion_tensor ( const cs_mesh_t * m,
int idiffp,
double thetap,
const cs_real_66_t cofbfp[],
const cs_real_66_t fimp[],
const cs_real_66_t i_visc[],
const cs_real_t b_visc[],
cs_real_66_t *restrict da,
cs_real_66_t *restrict xa )

Build the diffusion matrix for a tensor field with a tensorial diffusivity (symmetric matrix).

Build the diffusion matrix for a vector field with a tensorial diffusivity (symmetric matrix).

The diffusion is not reconstructed. The matrix is split into a diagonal block (3x3 times number of cells) and an extra diagonal part (of dimension 3x3 the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]cofbfpboundary condition array for the variable flux (implicit part - 3x3 tensor array)
[in]fimp$ \tens{f_s}^{imp} $
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_sym_matrix_scalar()

void cs_sym_matrix_scalar ( const cs_mesh_t * m,
int idiffp,
double thetap,
const cs_real_t cofbfp[],
const cs_real_t rovsdt[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_t *restrict da,
cs_real_t *restrict xa )

Build the diffusion matrix for a scalar field. (symmetric matrix).

The diffusion is not reconstructed. The matrix is split into a diagonal block (number of cells) and an extra diagonal part (of dimension the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]cofbfpboundary condition array for the variable flux (implicit part)
[in]rovsdtworking array
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ S_\fib $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra diagonal part of the matrix

◆ cs_sym_matrix_tensor()

void cs_sym_matrix_tensor ( const cs_mesh_t * m,
int idiffp,
double thetap,
const cs_real_66_t cofbfts[],
const cs_real_66_t fimp[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_66_t *restrict da,
cs_real_t *restrict xa )

Build the diffusion matrix for a tensor field (symmetric matrix).

The diffusion is not reconstructed. The matrix is split into a diagonal block (6x6 times number of cells) and an extra diagonal part (of dimension the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]cofbftsboundary condition array for the variable flux (Implicit part - 6x6 tensor array)
[in]fimppart of the diagonal
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ cs_sym_matrix_vector()

void cs_sym_matrix_vector ( const cs_mesh_t * m,
int idiffp,
double thetap,
const cs_real_33_t cofbfp[],
const cs_real_33_t fimp[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_33_t *restrict da,
cs_real_t *restrict xa )

Build the diffusion matrix for a vector field (symmetric matrix).

The diffusion is not reconstructed. The matrix is split into a diagonal block (3x3 times number of cells) and an extra diagonal part (of dimension the number of internal faces).

Parameters
[in]mpointer to mesh structure
[in]idiffpindicator
  • 1 diffusion
  • 0 otherwise
[in]thetapweighting coefficient for the theta-scheme,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]cofbfpboundary condition array for the variable flux (implicit part - 3x3 tensor array)
[in]fimp$ \tens{f_s}^{imp} $
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[out]dadiagonal part of the matrix
[out]xaextra interleaved diagonal part of the matrix

◆ matrdt()

void matrdt ( const int *const iconvp,
const int *const idiffp,
const int *const isym,
const cs_real_t coefbp[],
const cs_real_t cofbfp[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
cs_real_t da[] )

◆ matrix()

void matrix ( const int * iconvp,
const int * idiffp,
const int * ndircp,
const int * isym,
const cs_real_t * thetap,
const int * imucpp,
const cs_real_t coefbp[],
const cs_real_t cofbfp[],
const cs_real_t rovsdt[],
const cs_real_t i_massflux[],
const cs_real_t b_massflux[],
const cs_real_t i_visc[],
const cs_real_t b_visc[],
const cs_real_t xcpp[],
cs_real_t da[],
cs_real_t xa[] )