Skip to content

Commit

Permalink
start fixing tests
Browse files Browse the repository at this point in the history
  • Loading branch information
jdannberg committed Mar 27, 2018
1 parent 6b0d37d commit a1865a3
Show file tree
Hide file tree
Showing 85 changed files with 5,785 additions and 6,006 deletions.
8 changes: 5 additions & 3 deletions include/aspect/melt.h
Expand Up @@ -381,10 +381,12 @@ namespace aspect

/**
* Given the Darcy coefficient as computed by the material model, limit the
* coefficient to a minimum value based on the reference Darcy coefficient,
* and return this value.
* coefficient to a minimum value based on the reference Darcy coefficient
* in melt cells, and set it to zero in cells that are not melt cells, and
* return this value.
*/
double limited_darcy_coefficient(const double K_D) const;
double limited_darcy_coefficient(const double K_D,
const bool is_melt_cell) const;

/**
* The porosity limit for melt migration. For smaller porosities, the equations
Expand Down
2 changes: 1 addition & 1 deletion source/postprocess/visualization/melt.cc
Expand Up @@ -147,7 +147,7 @@ namespace aspect
}
else if (property_names[i] == "darcy coefficient")
{
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q]);
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q], p_c_scale > 0);
computed_quantities[q][output_index] = K_D;
}
else if (property_names[i] == "darcy coefficient no cutoff")
Expand Down
12 changes: 7 additions & 5 deletions source/simulator/melt.cc
Expand Up @@ -212,7 +212,7 @@ namespace aspect
const double eta = scratch.material_model_outputs.viscosities[q];
const double one_over_eta = 1. / eta;
const double eta_two_thirds = scratch.material_model_outputs.viscosities[q] * 2.0 / 3.0;
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q]);
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q], p_c_scale > 0);

/*
- R = 1/eta M_p + K_D L_p for p
Expand Down Expand Up @@ -433,7 +433,7 @@ namespace aspect
const unsigned int porosity_index = introspection.compositional_index_for_name("porosity");
const double porosity = std::max(scratch.material_model_inputs.composition[q][porosity_index],0.000);

const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q]);
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q], p_c_scale > 0);
const double viscosity_c = melt_outputs->compaction_viscosities[q];
const Tensor<1,dim> density_gradient_f = melt_outputs->fluid_density_gradients[q];
const double density_f = melt_outputs->fluid_densities[q];
Expand Down Expand Up @@ -535,6 +535,7 @@ namespace aspect
const FEValuesExtractors::Scalar ex_p_f = introspection.variable("fluid pressure").extractor_scalar();
const unsigned int p_f_component_index = introspection.variable("fluid pressure").first_component_index;
const unsigned int p_c_component_index = introspection.variable("compaction pressure").first_component_index;
const double p_c_scale = dynamic_cast<const MaterialModel::MeltInterface<dim>*>(&this->get_material_model())->p_c_scale(scratch.material_model_inputs, scratch.material_model_outputs, this->get_melt_handler(), true);

const typename DoFHandler<dim>::face_iterator face = scratch.face_material_model_inputs.current_cell->face(scratch.face_number);

Expand All @@ -557,7 +558,7 @@ namespace aspect
const Tensor<1,dim>
gravity = this->get_gravity_model().gravity_vector (scratch.face_finite_element_values.quadrature_point(q));
const double density_f = melt_outputs->fluid_densities[q];
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q]);
const double K_D = this->get_melt_handler().limited_darcy_coefficient(melt_outputs->permeabilities[q] / melt_outputs->fluid_viscosities[q], p_c_scale > 0);

for (unsigned int i=0, i_stokes=0; i_stokes<stokes_dofs_per_cell; /*increment at end of loop*/)
{
Expand Down Expand Up @@ -1590,10 +1591,11 @@ namespace aspect
template <int dim>
double
MeltHandler<dim>::
limited_darcy_coefficient(const double K_D) const
limited_darcy_coefficient(const double K_D,
const bool is_melt_cell) const
{
const double ref_K_D = dynamic_cast<const MaterialModel::MeltInterface<dim>*>(&this->get_material_model())->reference_darcy_coefficient();
return std::max(K_D, 1.e-3*ref_K_D);
return is_melt_cell ? std::max(K_D, 1.e-3*ref_K_D) : 0;
}


Expand Down
24 changes: 12 additions & 12 deletions tests/adiabatic_heating_with_melt/screen-output
@@ -1,13 +1,13 @@

Number of active cells: 320 (on 3 levels)
Number of degrees of freedom: 9,270 (2,898+810+2,898+405+1,449+405+405)
Number of degrees of freedom: 9,825 (2,898+1,365+2,898+405+1,449+405+405)

*** Timestep 0: t=0 seconds
Solving temperature system... 0 iterations.
Solving porosity system ... 0 iterations.
Skipping peridotite composition solve because RHS is zero.
Rebuilding Stokes preconditioner...
Solving Stokes system... 54+0 iterations.
Solving Stokes system... 58+0 iterations.
Solving for u_f in 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.85418e-16, 9.58741e-17, 0, 1
Relative nonlinear residual (total system) after nonlinear iteration 1: 1
Expand All @@ -18,8 +18,8 @@ Number of degrees of freedom: 9,270 (2,898+810+2,898+405+1,449+405+405)
Rebuilding Stokes preconditioner...
Solving Stokes system... 5+0 iterations.
Solving for u_f in 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.85418e-16, 9.58741e-17, 0, 9.29021e-13
Relative nonlinear residual (total system) after nonlinear iteration 2: 9.29021e-13
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.85418e-16, 9.58741e-17, 0, 3.92921e-13
Relative nonlinear residual (total system) after nonlinear iteration 2: 3.92921e-13


Postprocessing:
Expand All @@ -31,9 +31,9 @@ Number of degrees of freedom: 9,270 (2,898+810+2,898+405+1,449+405+405)
Solving porosity system ... 0 iterations.
Skipping peridotite composition solve because RHS is zero.
Rebuilding Stokes preconditioner...
Solving Stokes system... 17+0 iterations.
Solving Stokes system... 18+0 iterations.
Solving for u_f in 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.0391e-06, 1.53408e-16, 0, 4.87025e-09
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.0391e-06, 5.68755e-14, 0, 6.85769e-09
Relative nonlinear residual (total system) after nonlinear iteration 1: 2.0391e-06


Expand All @@ -43,12 +43,12 @@ Number of degrees of freedom: 9,270 (2,898+810+2,898+405+1,449+405+405)

*** Timestep 2: t=2.29425e+13 seconds
Solving temperature system... 13 iterations.
Solving porosity system ... 0 iterations.
Solving porosity system ... 2 iterations.
Skipping peridotite composition solve because RHS is zero.
Rebuilding Stokes preconditioner...
Solving Stokes system... 20+0 iterations.
Solving Stokes system... 22+0 iterations.
Solving for u_f in 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.21934e-06, 1.21242e-16, 0, 4.20135e-09
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.21934e-06, 6.65311e-10, 0, 5.93641e-09
Relative nonlinear residual (total system) after nonlinear iteration 1: 2.21934e-06


Expand All @@ -58,12 +58,12 @@ Number of degrees of freedom: 9,270 (2,898+810+2,898+405+1,449+405+405)

*** Timestep 3: t=2.5e+13 seconds
Solving temperature system... 6 iterations.
Solving porosity system ... 0 iterations.
Solving porosity system ... 1 iterations.
Skipping peridotite composition solve because RHS is zero.
Rebuilding Stokes preconditioner...
Solving Stokes system... 16+0 iterations.
Solving Stokes system... 18+0 iterations.
Solving for u_f in 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 4.7975e-08, 1.56341e-16, 0, 8.51551e-10
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 4.7975e-08, 3.90688e-11, 0, 1.19967e-09
Relative nonlinear residual (total system) after nonlinear iteration 1: 4.7975e-08


Expand Down
8 changes: 4 additions & 4 deletions tests/adiabatic_heating_with_melt/statistics
Expand Up @@ -20,7 +20,7 @@
# 20: Total adiabatic heating rate (W)
# 21: Average adiabatic heating of melt rate (W/kg)
# 22: Total adiabatic heating of melt rate (W)
0 0.000000000000e+00 0.000000000000e+00 320 3303 1449 810 2 0 0 0 59 62 62 1.72315000e+03 1.75516666e+03 1.78757747e+03 8.48573063e-03 -3.19121211e-11 -1.77378538e+02 3.19121211e-11 1.77378538e+02
1 1.147124999994e+13 1.147124999994e+13 320 3303 1449 810 1 8 0 0 17 18 18 1.72319196e+03 1.75516670e+03 1.78757747e+03 8.47462138e-03 -3.19121219e-11 -1.77378542e+02 -1.02517773e-26 -5.69829019e-14
2 2.294249984878e+13 1.147124984883e+13 320 3303 1449 810 1 13 0 0 20 21 21 1.72316615e+03 1.75516665e+03 1.78757747e+03 8.48144718e-03 -3.19121210e-11 -1.77378537e+02 2.35199516e-27 1.30731975e-14
3 2.500000000000e+13 2.057500151224e+12 320 3303 1449 810 1 6 0 0 16 17 17 1.72315936e+03 1.75516664e+03 1.78757747e+03 8.48324430e-03 -3.19121208e-11 -1.77378536e+02 4.13359881e-26 2.29759630e-13
0 0.000000000000e+00 0.000000000000e+00 320 3303 1449 810 2 0 0 0 63 66 461 1.72315000e+03 1.75516666e+03 1.78757747e+03 8.48573063e-03 -3.19121211e-11 -1.77378538e+02 3.19121211e-11 1.77378538e+02
1 1.147124999996e+13 1.147124999996e+13 320 3303 1449 810 1 8 0 0 18 19 133 1.72319196e+03 1.75516670e+03 1.78757747e+03 8.47462138e-03 -3.19121219e-11 -1.77378542e+02 -9.46965197e-27 -5.26355805e-14
2 2.294249984879e+13 1.147124984883e+13 320 3303 1449 810 1 13 2 0 22 23 160 1.72316615e+03 1.75516665e+03 1.78757747e+03 8.48144718e-03 -3.19121210e-11 -1.77378537e+02 2.48973431e-27 1.38387990e-14
3 2.500000000000e+13 2.057500151206e+12 320 3303 1449 810 1 6 1 0 18 19 132 1.72315936e+03 1.75516664e+03 1.78757747e+03 8.48324430e-03 -3.19121208e-11 -1.77378536e+02 4.93577723e-26 2.74347464e-13
14 changes: 6 additions & 8 deletions tests/compaction_length_refinement/screen-output
@@ -1,8 +1,6 @@
-----------------------------------------------------------------------------
-----------------------------------------------------------------------------

Number of active cells: 1,024 (on 6 levels)
Number of degrees of freedom: 32,842 (8,450+2,178+8,450+1,089+4,225+4,225+4,225)
Number of degrees of freedom: 34,825 (8,450+4,161+8,450+1,089+4,225+4,225+4,225)

*** Timestep 0: t=0 years
Solving temperature system... 0 iterations.
Expand All @@ -16,7 +14,7 @@ Number of degrees of freedom: 32,842 (8,450+2,178+8,450+1,089+4,225+4,225+4,225)


Number of active cells: 1,117 (on 7 levels)
Number of degrees of freedom: 36,960 (9,510+2,450+9,510+1,225+4,755+4,755+4,755)
Number of degrees of freedom: 39,086 (9,510+4,576+9,510+1,225+4,755+4,755+4,755)

*** Timestep 0: t=0 years
Solving temperature system... 0 iterations.
Expand All @@ -30,7 +28,7 @@ Number of degrees of freedom: 36,960 (9,510+2,450+9,510+1,225+4,755+4,755+4,755)


Number of active cells: 1,696 (on 8 levels)
Number of degrees of freedom: 56,069 (14,434+3,700+14,434+1,850+7,217+7,217+7,217)
Number of degrees of freedom: 59,307 (14,434+6,938+14,434+1,850+7,217+7,217+7,217)

*** Timestep 0: t=0 years
Solving temperature system... 0 iterations.
Expand All @@ -44,7 +42,7 @@ Number of degrees of freedom: 56,069 (14,434+3,700+14,434+1,850+7,217+7,217+7,21


Number of active cells: 2,023 (on 9 levels)
Number of degrees of freedom: 67,005 (17,250+4,420+17,250+2,210+8,625+8,625+8,625)
Number of degrees of freedom: 70,864 (17,250+8,279+17,250+2,210+8,625+8,625+8,625)

*** Timestep 0: t=0 years
Solving temperature system... 0 iterations.
Expand All @@ -58,7 +56,7 @@ Number of degrees of freedom: 67,005 (17,250+4,420+17,250+2,210+8,625+8,625+8,62


Number of active cells: 2,026 (on 9 levels)
Number of degrees of freedom: 67,091 (17,272+4,426+17,272+2,213+8,636+8,636+8,636)
Number of degrees of freedom: 70,956 (17,272+8,291+17,272+2,213+8,636+8,636+8,636)

*** Timestep 0: t=0 years
Solving temperature system... 0 iterations.
Expand All @@ -72,7 +70,7 @@ Number of degrees of freedom: 67,091 (17,272+4,426+17,272+2,213+8,636+8,636+8,63


Number of active cells: 2,023 (on 9 levels)
Number of degrees of freedom: 67,005 (17,250+4,420+17,250+2,210+8,625+8,625+8,625)
Number of degrees of freedom: 70,864 (17,250+8,279+17,250+2,210+8,625+8,625+8,625)

*** Timestep 0: t=0 years
Solving temperature system... 0 iterations.
Expand Down
1 change: 0 additions & 1 deletion tests/composition_active_with_melt.cc
Expand Up @@ -75,7 +75,6 @@ namespace aspect

for (unsigned int i=0; i<in.position.size(); ++i)
{
double porosity = in.composition[i][0];
melt_out->compaction_viscosities[i] = 1.0;
melt_out->fluid_viscosities[i] = 1.0;
melt_out->permeabilities[i] = 0.0; //1e-30*porosity * porosity;
Expand Down
27 changes: 6 additions & 21 deletions tests/composition_active_with_melt.prm
Expand Up @@ -10,9 +10,11 @@ set Dimension = 2
set Start time = 0
set End time = 1
set Use years in output instead of seconds = false
set Number of cheap Stokes solver steps = 0
set Linear solver tolerance = 1e-8
#set Nonlinear solver scheme = iterated IMPES

subsection Discretization
set Stokes velocity polynomial degree = 2
set Composition polynomial degree = 1
end

subsection Geometry model
set Model name = box
Expand All @@ -26,7 +28,6 @@ end

subsection Model settings
set Fixed temperature boundary indicators = 2, 3
# set Fixed composition boundary indicators = 2, 3
set Zero velocity boundary indicators =
set Tangential velocity boundary indicators = 0, 1, 2
set Prescribed velocity boundary indicators = 3: function
Expand All @@ -36,6 +37,7 @@ end

subsection Melt settings
set Melt transport threshold = 0.0
set Use discontinuous compaction pressure = false
end

subsection Heating model
Expand Down Expand Up @@ -115,20 +117,3 @@ subsection Initial composition model
set Function expression = if(y<0.2, 1, 0) ; 0.0
end
end

subsection Discretization
set Stokes velocity polynomial degree = 2
set Composition polynomial degree = 1
end

subsection Boundary fluid pressure model
set Plugin name = density
subsection Density
set Density formulation = solid density
end

end

subsection Boundary composition model
set List of model names = initial composition
end

0 comments on commit a1865a3

Please sign in to comment.