Commit 70923058 authored by Jonathan Lambrechts's avatar Jonathan Lambrechts
Browse files

re-use a_cg

parent 1303ba9c
Pipeline #4929 passed with stage
in 26 seconds
......@@ -695,27 +695,11 @@ static void fluid_problem_volume(FluidProblem *problem, const double *solution_o
free(h1);
}
static void fluid_problem_surface_tension_bnd(FluidProblem *problem, double *all_local_vector){
static void fluid_problem_surface_tension_bnd(FluidProblem *problem, double *a_cg, double *all_local_vector){
const Mesh *mesh = problem->mesh;
double sigma = problem->sigma;
const int n_fields = fluid_problem_n_fields(problem);
const size_t local_size = N_SF*n_fields;
double *a_cg = malloc(sizeof(double)*mesh->n_nodes);
for (int i = 0; i< mesh->n_nodes; ++i) {
a_cg[i] = 0;
}
for (int iel = 0; iel < mesh->n_elements; ++iel) {
const int *el = &mesh->elements[iel*N_N];
double dxidx[D][D];
const double det = mesh_dxidx(mesh,iel,dxidx);
const double vol = node_volume_from_detj(det);
for (int i = 0; i< N_N; ++i) {
a_cg[el[i]] += vol * problem->concentration[iel*N_N+i];
}
}
for (int i = 0; i< mesh->n_nodes; ++i) {
a_cg[i] /= problem->node_volume[i];
}
for (int ibnd = 0; ibnd < problem->n_weak_boundaries; ++ibnd){
WeakBoundary *wbnd = problem->weak_boundaries + ibnd;
int bndid= -1;
......@@ -777,30 +761,13 @@ static void fluid_problem_surface_tension_bnd(FluidProblem *problem, double *all
}
}
}
free(a_cg);
}
static void fluid_problem_surface_tension(FluidProblem *problem, double *all_local_vector){
static void fluid_problem_surface_tension(FluidProblem *problem, double *a_cg, double *all_local_vector){
if (problem->n_fluids == 1) return;
const Mesh *mesh = problem->mesh;
double sigma = problem->sigma;
int n_fields = fluid_problem_n_fields(problem);
double *a_cg = malloc(sizeof(double)*mesh->n_nodes);
for (int i = 0; i< mesh->n_nodes; ++i) {
a_cg[i] = 0;
}
for (int iel = 0; iel < mesh->n_elements; ++iel) {
const int *el = &mesh->elements[iel*N_N];
double dxidx[D][D];
const double det = mesh_dxidx(mesh,iel,dxidx);
const double vol = node_volume_from_detj(det);
for (int i = 0; i< N_N; ++i) {
a_cg[el[i]] += vol * problem->concentration[iel*N_N+i];
}
}
for (int i = 0; i< mesh->n_nodes; ++i) {
a_cg[i] /= problem->node_volume[i];
}
for (int iel = 0; iel < mesh->n_elements; ++iel) {
const int *el = &mesh->elements[iel*N_N];
double dxidx[D][D], dphi[N_N][D];
......@@ -830,7 +797,26 @@ static void fluid_problem_surface_tension(FluidProblem *problem, double *all_loc
}
}
}
free(a_cg);
}
static void fluid_problem_dg_to_cg(FluidProblem *problem, const double *dg, double *cg) {
const Mesh *mesh = problem->mesh;
for (int i = 0; i< mesh->n_nodes; ++i) {
cg[i] = 0;
}
for (int iel = 0; iel < mesh->n_elements; ++iel) {
const int *el = &mesh->elements[iel*N_N];
double dxidx[D][D];
const double det = mesh_dxidx(mesh,iel,dxidx);
const double vol = node_volume_from_detj(det);
for (int i = 0; i< N_N; ++i) {
cg[el[i]] += vol * dg[iel*N_N+i];
}
}
for (int i = 0; i< mesh->n_nodes; ++i) {
cg[i] /= problem->node_volume[i];
}
}
static void fluid_problem_assemble_system(FluidProblem *problem, double *rhs, const double *solution_old, double dt)
......@@ -871,8 +857,11 @@ static void fluid_problem_assemble_system(FluidProblem *problem, double *rhs, co
compute_weak_boundary_conditions(problem, dt, all_local_vector, all_local_matrix);
fluid_problem_volume(problem, solution_old, dt, all_local_vector, all_local_matrix);
if (problem->n_fluids == 2){
fluid_problem_surface_tension(problem, all_local_vector);
fluid_problem_surface_tension_bnd(problem, all_local_vector);
double *a_cg = malloc(sizeof(double)*mesh->n_nodes);
fluid_problem_dg_to_cg(problem, problem->concentration, a_cg);
fluid_problem_surface_tension(problem, a_cg, all_local_vector);
fluid_problem_surface_tension_bnd(problem, a_cg, all_local_vector);
free(a_cg);
}
for (int iel=0; iel < mesh->n_elements; ++iel) {
double *local_vector = &all_local_vector[local_size*iel];
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment