#include "cs_defs.h"
#include <chrono>
#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_dispatch.h"
#include "cs_halo.h"
#include "cs_halo_perio.h"
#include "cs_log.h"
#include "cs_matrix.h"
#include "cs_matrix_default.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"
Functions | |
void | cs_matrix_compute_coeffs (cs_matrix_t *a, const cs_field_t *f, int iconvp, int idiffp, int ndircp, double thetap, double relaxp, int imucpp, const cs_field_bc_coeffs_t *bc_coeffs, 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[]) |
Build the diagonal of the advection/diffusion matrix for determining the variable time step, flow, Fourier. More... | |
template<cs_lnum_t stride> | |
void | cs_matrix_compute_coeffs (cs_matrix_t *a, const cs_field_t *f, int iconvp, int idiffp, int tensorial_diffusion, int ndircp, cs_lnum_t eb_size, double thetap, double relaxp, const cs_field_bc_coeffs_t *bc_coeffs, const cs_real_t fimp[][stride][stride], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[]) |
Build the diagonal of the advection/diffusion matrix for determining the variable time step, flow, Fourier. More... | |
template void | cs_matrix_compute_coeffs (cs_matrix_t *a, const cs_field_t *f, int iconvp, int idiffp, int tensorial_diffusion, int ndircp, cs_lnum_t eb_size, double thetap, double relaxp, const cs_field_bc_coeffs_t *bc_coeffs, const cs_real_t fimp[][3][3], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[]) |
template void | cs_matrix_compute_coeffs (cs_matrix_t *a, const cs_field_t *f, int iconvp, int idiffp, int tensorial_diffusion, int ndircp, cs_lnum_t eb_size, double thetap, double relaxp, const cs_field_bc_coeffs_t *bc_coeffs, const cs_real_t fimp[][6][6], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[]) |
void | cs_matrix_wrapper (int iconvp, int idiffp, int ndircp, int isym, double thetap, int imucpp, const cs_field_bc_coeffs_t *bc_coeffs, 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[]) |
template<cs_lnum_t stride> | |
void | cs_matrix_wrapper (int iconvp, int idiffp, int tensorial_diffusion, int ndircp, int isym, cs_lnum_t eb_size, double thetap, const cs_field_bc_coeffs_t *bc_coeffs_v, const cs_real_t fimp[][stride][stride], 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[][stride][stride], cs_real_t xa[]) |
Build the matrix for a vector or tensor field. More... | |
template void | cs_matrix_wrapper (int iconvp, int idiffp, int tensorial_diffusion, int ndircp, int isym, cs_lnum_t eb_size, double thetap, const cs_field_bc_coeffs_t *bc_coeffs_v, const cs_real_t fimp[][3][3], 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[][3][3], cs_real_t xa[]) |
template void | cs_matrix_wrapper (int iconvp, int idiffp, int tensorial_diffusion, int ndircp, int isym, cs_lnum_t eb_size, double thetap, const cs_field_bc_coeffs_t *bc_coeffs_v, const cs_real_t fimp[][6][6], 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[][6][6], cs_real_t xa[]) |
void | cs_matrix_time_step (const cs_mesh_t *m, int iconvp, int idiffp, int isym, const cs_field_bc_coeffs_t *bc_coeffs, 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. More... | |
void cs_matrix_compute_coeffs | ( | cs_matrix_t * | a, |
const cs_field_t * | f, | ||
int | iconvp, | ||
int | idiffp, | ||
int | ndircp, | ||
double | thetap, | ||
double | relaxp, | ||
int | imucpp, | ||
const cs_field_bc_coeffs_t * | bc_coeffs, | ||
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[] | ||
) |
Build the diagonal of the advection/diffusion matrix for determining the variable time step, flow, Fourier.
[in,out] | a | pointer to matrix structure |
[in] | f | pointer to field, or null |
[in] | iconvp | indicator
|
[in] | idiffp | indicator
|
[in] | ndircp | number of Dirichlet BCs |
[in] | thetap | time scheme parameter |
[in] | relaxp | relaxation coefficient (if < 1) |
[in] | imucp | 1 for temperature (with Cp), 0 otherwise |
[in] | bc_coeffs | boundary condition structure |
[in] | rovsdt | implicit terms (rho / dt) |
[in] | i_massflux | mass flux at interior faces |
[in] | b_massflux | mass flux at border faces |
[in] | i_visc | ![]() |
[in] | b_visc | ![]() |
[in] | xcpp | Cp per cell, or null |
template void cs_matrix_compute_coeffs | ( | cs_matrix_t * | a, |
const cs_field_t * | f, | ||
int | iconvp, | ||
int | idiffp, | ||
int | tensorial_diffusion, | ||
int | ndircp, | ||
cs_lnum_t | eb_size, | ||
double | thetap, | ||
double | relaxp, | ||
const cs_field_bc_coeffs_t * | bc_coeffs, | ||
const cs_real_t | fimp[][3][3], | ||
const cs_real_t | i_massflux[], | ||
const cs_real_t | b_massflux[], | ||
const cs_real_t | i_visc[], | ||
const cs_real_t | b_visc[] | ||
) |
template void cs_matrix_compute_coeffs | ( | cs_matrix_t * | a, |
const cs_field_t * | f, | ||
int | iconvp, | ||
int | idiffp, | ||
int | tensorial_diffusion, | ||
int | ndircp, | ||
cs_lnum_t | eb_size, | ||
double | thetap, | ||
double | relaxp, | ||
const cs_field_bc_coeffs_t * | bc_coeffs, | ||
const cs_real_t | fimp[][6][6], | ||
const cs_real_t | i_massflux[], | ||
const cs_real_t | b_massflux[], | ||
const cs_real_t | i_visc[], | ||
const cs_real_t | b_visc[] | ||
) |
void cs_matrix_compute_coeffs | ( | cs_matrix_t * | a, |
const cs_field_t * | f, | ||
int | iconvp, | ||
int | idiffp, | ||
int | tensorial_diffusion, | ||
int | ndircp, | ||
cs_lnum_t | eb_size, | ||
double | thetap, | ||
double | relaxp, | ||
const cs_field_bc_coeffs_t * | bc_coeffs, | ||
const cs_real_t | fimp[][stride][stride], | ||
const cs_real_t | i_massflux[], | ||
const cs_real_t | b_massflux[], | ||
const cs_real_t | i_visc[], | ||
const cs_real_t | b_visc[] | ||
) |
Build the diagonal of the advection/diffusion matrix for determining the variable time step, flow, Fourier.
stride | 3 for vectors, 6 for tensors |
[in,out] | a | pointer to matrix structure |
[in] | f | pointer to field, or null |
[in] | iconvp | indicator
|
[in] | idiffp | indicator
|
[in] | tensorial_diffusion | indicator |
[in] | ndircp | number of Dirichlet BCs |
[in] | thetap | time scheme parameter |
[in] | relaxp | relaxation coefficient (if < 1) |
[in] | eb_size | extra-diagonal block size (1 or 3 for stride 3, 1 for stride 6) |
[in] | bc_coeffs | boundary conditions structure |
[in] | fimp | implicit terms, or null |
[in] | i_massflux | mass flux at interior faces |
[in] | b_massflux | mass flux at border faces |
[in] | i_visc | ![]() |
[in] | b_visc | ![]() |
void cs_matrix_time_step | ( | const cs_mesh_t * | m, |
int | iconvp, | ||
int | idiffp, | ||
int | isym, | ||
const cs_field_bc_coeffs_t * | bc_coeffs, | ||
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.
[in] | m | pointer to mesh structure |
[in] | iconvp | indicator
|
[in] | idiffp | indicator
|
[in] | isym | indicator
|
[in] | bc_coeffs | boundary condition structure for the variable |
[in] | i_massflux | mass flux at interior faces |
[in] | b_massflux | mass flux at border faces |
[in] | i_visc | ![]() |
[in] | b_visc | ![]() |
[out] | da | diagonal part of the matrix |
void cs_matrix_wrapper | ( | int | iconvp, |
int | idiffp, | ||
int | ndircp, | ||
int | isym, | ||
double | thetap, | ||
int | imucpp, | ||
const cs_field_bc_coeffs_t * | bc_coeffs, | ||
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[] | ||
) |
template void cs_matrix_wrapper | ( | int | iconvp, |
int | idiffp, | ||
int | tensorial_diffusion, | ||
int | ndircp, | ||
int | isym, | ||
cs_lnum_t | eb_size, | ||
double | thetap, | ||
const cs_field_bc_coeffs_t * | bc_coeffs_v, | ||
const cs_real_t | fimp[][3][3], | ||
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[][3][3], | ||
cs_real_t | xa[] | ||
) |
template void cs_matrix_wrapper | ( | int | iconvp, |
int | idiffp, | ||
int | tensorial_diffusion, | ||
int | ndircp, | ||
int | isym, | ||
cs_lnum_t | eb_size, | ||
double | thetap, | ||
const cs_field_bc_coeffs_t * | bc_coeffs_v, | ||
const cs_real_t | fimp[][6][6], | ||
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[][6][6], | ||
cs_real_t | xa[] | ||
) |
void cs_matrix_wrapper | ( | int | iconvp, |
int | idiffp, | ||
int | tensorial_diffusion, | ||
int | ndircp, | ||
int | isym, | ||
cs_lnum_t | eb_size, | ||
double | thetap, | ||
const cs_field_bc_coeffs_t * | bc_coeffs_v, | ||
const cs_real_t | fimp[][stride][stride], | ||
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[][stride][stride], | ||
cs_real_t | xa[] | ||
) |
Build the matrix for a vector or tensor field.
The advection (if present) is upwind. The diffusion is not reconstructed. The matrix is split into a diagonal part (stride*stride blocks) and an extra diagonal part.
stride | 3 for vectors, 6 for tensors |
[in] | m | pointer to mesh structure |
[in] | idiffp | indicator
|
[in] | thetap | weighting coefficient for the theta-scheme,
|
[in] | bc_coeffs_v | boundary condition structure for the variable |
[in] | fimp | ![]() |
[in] | i_visc | ![]() |
[in] | b_visc | ![]() |
[out] | da | diagonal part of the matrix |
[out] | xa | extra interleaved diagonal part of the matrix |