diff --git a/examples/gringarten_multiapp/frect1.i b/examples/gringarten_multiapp/frect1.i index 8e8680dd..855c0ceb 100644 --- a/examples/gringarten_multiapp/frect1.i +++ b/examples/gringarten_multiapp/frect1.i @@ -1,17 +1,29 @@ # Cold water injection into one side of the fracture network, and production from the other side -injection_rate = 0.1 #25 # kg/s -endTime = 3.16e8 +injection_rate = 0.1 # kg/s +endTime = 3e8 [Mesh] uniform_refine = 0 - [single_frac] - type = FileMeshGenerator - file = 'rect_100_10_medium.e' + [generate] + type = GeneratedMeshGenerator + dim = 2 + nx = 10 + xmin = -5 + xmax = 5 + ny = 50 + ymin = -50 + ymax = 50 + [] + [rotate] + type = TransformGenerator + input = generate + transform = ROTATE + vector_value = '90 90 90' [] [injection_node] type = BoundingBoxNodeSetGenerator - input = single_frac - bottom_left = '-5 -55 -5' - top_right = '5 -45 5' + input = rotate + bottom_left = '-0.05 -51 -5.1' + top_right = '0.05 -49.9 5.1' new_boundary = injection_node [] [] @@ -45,8 +57,9 @@ endTime = 3.16e8 coupling_type = ThermoHydro porepressure = frac_P temperature = frac_T - fp = water + fp = the_simple_fluid #water pressure_unit = Pa + # stabilization = full [] [Kernels] @@ -129,10 +142,10 @@ endTime = 3.16e8 [heat_transfer_coefficient_auxk] type = ParsedAux variable = heat_transfer_coefficient - args = 'enclosing_element_normal_length enclosing_element_normal_thermal_cond' + coupled_variables = 'enclosing_element_normal_length enclosing_element_normal_thermal_cond' constant_names = h_s constant_expressions = 1E3 #This is the value being assigned to h_s. Should be much bigger than thermal_conductivity / L ~ 1 - function = 'if(enclosing_element_normal_length = 0, 0, h_s * enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length / (h_s * enclosing_element_normal_length * enclosing_element_normal_length + enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length))' + expression = 'if(enclosing_element_normal_length = 0, 0, h_s * enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length / (h_s * enclosing_element_normal_length * enclosing_element_normal_length + enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length))' [] [insitu_pp] type = FunctionAux @@ -219,11 +232,17 @@ endTime = 3.16e8 [] [FluidProperties] + [the_simple_fluid] + type = SimpleFluidProperties + bulk_modulus = 2E9 + viscosity = 1.0E-3 + density0 = 1000.0 + [] [true_water] type = Water97FluidProperties [] [water] - type = TabulatedFluidProperties + type = TabulatedBicubicFluidProperties fp = true_water temperature_min = 275 # K temperature_max = 600 @@ -233,27 +252,21 @@ endTime = 3.16e8 [] [Materials] - [porosity] - type = PorousFlowPorosityLinear - porosity_ref = 1E-4 - P_ref = insitu_pp - P_coeff = 3e-10 - porosity_min = 1E-5 - [] - [permeability] - type = PorousFlowPermeabilityKozenyCarman - k0 = 1E-15 - poroperm_function = kozeny_carman_phi0 - m = 0 - n = 3 - phi0 = 1E-4 - [] - [internal_energy] + [porosity_frac] + type = PorousFlowPorosityConst + porosity = 0.01 + [] + [permeability_frac] + type = PorousFlowPermeabilityConst + permeability = '1e-12 0 0 0 1e-12 0 0 0 1e-12' + [] + + [internal_energy_frac] type = PorousFlowMatrixInternalEnergy density = 2700 specific_heat_capacity = 0 [] - [aq_thermal_conductivity] + [aq_thermal_conductivity_frac] type = PorousFlowThermalConductivityIdeal dry_thermal_conductivity = '0.6E-4 0 0 0 0.6E-4 0 0 0 0.6E-4' [] @@ -262,17 +275,17 @@ endTime = 3.16e8 [Functions] [kg_rate] type = ParsedFunction - vals = 'dt kg_out' - vars = 'dt kg_out' - value = 'kg_out/dt' + symbol_values = 'dt kg_out' + symbol_names = 'dt kg_out' + expression = 'kg_out/dt' [] [insitu_pp] type = ParsedFunction - value = '9.81*1000*(3000 - z)' + expression = '9.81*1000*(3000 - z)' [] [insitu_T] type = ParsedFunction - value = '363' + expression = '363' [] [] @@ -284,19 +297,23 @@ endTime = 3.16e8 [kg_out] type = PorousFlowPlotQuantity uo = kg_out_uo + outputs = none [] [kg_per_s] type = FunctionValuePostprocessor function = kg_rate + outputs = none [] [J_out] type = PorousFlowPlotQuantity uo = J_out_uo + outputs = none [] [TK_in] type = PointValue variable = frac_T point = '0 -50 0' + outputs = none [] [TK_out] type = PointValue @@ -307,36 +324,46 @@ endTime = 3.16e8 type = PointValue variable = frac_P point = '0 50 0' + outputs = none [] [P_in] type = PointValue variable = frac_P point = '0 -50 0' + outputs = none [] [] [VectorPostprocessors] [heat_transfer_rate] type = NodalValueSampler - outputs = none sort_by = id variable = joules_per_s + outputs = none [] [] [Preconditioning] - [./superlu] + [superlu] type = SMP full = true petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' petsc_options_value = 'gmres lu superlu_dist' - [../] + [] [] +[Functions] + [dts] + type = ParsedFunction + expression = if(t<1e6,t*t,1e6) + [] +[] [Executioner] type = Transient solve_type = NEWTON - dt = 0.1e7 + # dt = 1e6 + dtmin = 1 + dtmax = 1e6 end_time = ${endTime} line_search = 'none' automatic_scaling = true @@ -346,21 +373,14 @@ endTime = 3.16e8 nl_max_its = 20 nl_rel_tol = 5e-04 nl_abs_tol = 1e-09 + [TimeStepper] + type = FunctionDT + function = dts + min_dt = 2 + [] [] [Outputs] print_linear_residuals = false - exodus = false csv = true - # [fracCSV] - # type = CSV - # sync_times = '100 200 300 400 500 600 700 800 900 - # 1000 2000 3000 4000 5000 6000 7000 8000 9000 - # 1000e1 2000e1 3000e1 4000e1 5000e1 6000e1 7000e1 8000e1 9000e1 - # 1000e2 2000e2 3000e2 4000e2 5000e2 6000e2 7000e2 8000e2 9000e2 - # 1000e3 2000e3 3000e3 4000e3 5000e3 6000e3 7000e3 8000e3 9000e3 - # 1000e4 2000e4 3000e4 4000e4 5000e4 6000e4 7000e4 8000e4 9000e4 - # 1000e5 2000e5 3000e5 4000e5 5000e5 6000e5 7000e5 8000e5 9000e5' - # sync_only = true - # [] [] diff --git a/examples/gringarten_multiapp/frect1_amr.i b/examples/gringarten_multiapp/frect1_amr.i deleted file mode 100644 index b5ea5cae..00000000 --- a/examples/gringarten_multiapp/frect1_amr.i +++ /dev/null @@ -1,386 +0,0 @@ -# Cold water injection into one side of the fracture network, and production from the other side -injection_rate = 0.1 #25 # kg/s -endTime = 3.16e8 - - -[Mesh] - uniform_refine = 0 - [generate] - type = GeneratedMeshGenerator - dim = 2 - nx = 10 - xmin = -5 - xmax = 5 - ny = 50 - ymin = -50 - ymax = 50 - [] - [./rotate] - type = TransformGenerator - input = generate - transform = ROTATE - vector_value = '90 90 90' - [] - [./offset] - type = TransformGenerator - input = rotate - transform = TRANSLATE - vector_value = '0.1 0 0' - [] - [injection_node] - type = BoundingBoxNodeSetGenerator - input = offset - bottom_left = '0.05 -51 -5.1' - top_right = '0.15 -49.9 5.1' - new_boundary = injection_node - [] -[] - -[GlobalParams] - PorousFlowDictator = dictator - gravity = '0 0 -9.81' -[] - -[Variables] - [frac_P] - [] - [frac_T] - [] -[] - -[ICs] - [frac_P] - type = FunctionIC - variable = frac_P - function = insitu_pp - [] - [frac_T] - type = FunctionIC - variable = frac_T - function = insitu_T - [] -[] - -[PorousFlowFullySaturated] - coupling_type = ThermoHydro - porepressure = frac_P - temperature = frac_T - fp = water - pressure_unit = Pa -[] - -[Kernels] - [toMatrix] - type = PorousFlowHeatMassTransfer - variable = frac_T - v = transferred_matrix_T - transfer_coefficient = heat_transfer_coefficient - save_in = joules_per_s - [] -[] - -[AuxVariables] - [transferred_matrix_T] - initial_condition = 363 - [] - [heat_transfer_coefficient] - family = MONOMIAL - order = CONSTANT - initial_condition = 0.0 - [] - [joules_per_s] - [] - [aperture] - family = MONOMIAL - order = CONSTANT - [] - [perm_times_app] - family = MONOMIAL - order = CONSTANT - [] - [density] - family = MONOMIAL - order = CONSTANT - [] - [viscosity] - family = MONOMIAL - order = CONSTANT - [] - [insitu_pp] - [] - [normal_dirn_x] - family = MONOMIAL - order = CONSTANT - [] - [normal_dirn_y] - family = MONOMIAL - order = CONSTANT - [] - [normal_dirn_z] - family = MONOMIAL - order = CONSTANT - [] - [enclosing_element_normal_length] - family = MONOMIAL - order = CONSTANT - [] - [enclosing_element_normal_thermal_cond] - family = MONOMIAL - order = CONSTANT - [] -[] - -[AuxKernels] - [normal_dirn_x_auxk] - type = PorousFlowElementNormal - variable = normal_dirn_x - component = x - [] - [normal_dirn_y] - type = PorousFlowElementNormal - variable = normal_dirn_y - component = y - [] - [normal_dirn_z] - type = PorousFlowElementNormal - variable = normal_dirn_z - component = z - [] - [heat_transfer_coefficient_auxk] - type = ParsedAux - variable = heat_transfer_coefficient - args = 'enclosing_element_normal_length enclosing_element_normal_thermal_cond' - constant_names = h_s - constant_expressions = 1E3 #This is the value being assigned to h_s. Should be much bigger than thermal_conductivity / L ~ 1 - function = 'if(enclosing_element_normal_length = 0, 0, h_s * enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length / (h_s * enclosing_element_normal_length * enclosing_element_normal_length + enclosing_element_normal_thermal_cond * 2 * enclosing_element_normal_length))' - [] - [insitu_pp] - type = FunctionAux - execute_on = initial - variable = insitu_pp - function = insitu_pp - [] - [aperture] - type = PorousFlowPropertyAux - variable = aperture - property = porosity - [] - [perm_times_app] - type = PorousFlowPropertyAux - variable = perm_times_app - property = permeability - row = 0 - column = 0 - [] - [density] - type = PorousFlowPropertyAux - variable = density - property = density - phase = 0 - [] - [viscosity] - type = PorousFlowPropertyAux - variable = viscosity - property = viscosity - phase = 0 - [] -[] - -[BCs] - [inject_heat] - type = DirichletBC - boundary = injection_node - variable = frac_T - value = 303 - [] -[] - -[DiracKernels] - [inject_fluid] - type = PorousFlowPointSourceFromPostprocessor - mass_flux = ${injection_rate} - point = '0.1 -50 0' - variable = frac_P - [] - [withdraw_fluid] - type = PorousFlowPeacemanBorehole - SumQuantityUO = kg_out_uo - bottom_p_or_t = 30.6e6 - character = 1 - line_length = 1 - point_file = production_single_fracture_amr.xyz - unit_weight = '0 0 0' - fluid_phase = 0 - use_mobility = true - variable = frac_P - [] - [withdraw_heat] - type = PorousFlowPeacemanBorehole - SumQuantityUO = J_out_uo - bottom_p_or_t = 30.6e6 - character = 1 - line_length = 1 - point_file = production_single_fracture_amr.xyz - unit_weight = '0 0 0' - fluid_phase = 0 - use_mobility = true - use_enthalpy = true - variable = frac_T - [] -[] - -[UserObjects] - [kg_out_uo] - type = PorousFlowSumQuantity - [] - [J_out_uo] - type = PorousFlowSumQuantity - [] -[] - -[FluidProperties] - [true_water] - type = Water97FluidProperties - [] - [water] - type = TabulatedFluidProperties - fp = true_water - temperature_min = 275 # K - temperature_max = 600 - interpolated_properties = 'density viscosity enthalpy internal_energy' - fluid_property_file = water97_tabulated.csv - [] -[] - -[Materials] - [porosity] - type = PorousFlowPorosityLinear - porosity_ref = 1E-4 - P_ref = insitu_pp - P_coeff = 3e-10 - porosity_min = 1E-5 - [] - [permeability] - type = PorousFlowPermeabilityKozenyCarman - k0 = 1E-15 - poroperm_function = kozeny_carman_phi0 - m = 0 - n = 3 - phi0 = 1E-4 - [] - [internal_energy] - type = PorousFlowMatrixInternalEnergy - density = 2700 - specific_heat_capacity = 0 - [] - [aq_thermal_conductivity] - type = PorousFlowThermalConductivityIdeal - dry_thermal_conductivity = '0.6E-4 0 0 0 0.6E-4 0 0 0 0.6E-4' - [] -[] - -[Functions] - [kg_rate] - type = ParsedFunction - vals = 'dt kg_out' - vars = 'dt kg_out' - value = 'kg_out/dt' - [] - [insitu_pp] - type = ParsedFunction - value = '9.81*1000*(3000 - z)' - [] - [insitu_T] - type = ParsedFunction - value = '363' - [] -[] - -[Postprocessors] - [dt] - type = TimestepSize - outputs = 'none' - [] - [kg_out] - type = PorousFlowPlotQuantity - uo = kg_out_uo - [] - [kg_per_s] - type = FunctionValuePostprocessor - function = kg_rate - [] - [J_out] - type = PorousFlowPlotQuantity - uo = J_out_uo - [] - [TK_in] - type = PointValue - variable = frac_T - point = '0.1 -50 0' - [] - [TK_out] - type = PointValue - variable = frac_T - point = '0.1 50 0' - [] - [P_in] - type = PointValue - variable = frac_P - point = '0.1 -50 0' - [] - [P_out] - type = PointValue - variable = frac_P - point = '0.1 50 0' - [] -[] - -[VectorPostprocessors] - [heat_transfer_rate] - type = NodalValueSampler - outputs = none - sort_by = id - variable = joules_per_s - [] -[] - -[Preconditioning] - [./superlu] - type = SMP - full = true - petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' - petsc_options_value = 'gmres lu superlu_dist' - [../] -[] - -[Executioner] - type = Transient - solve_type = NEWTON - dt = 0.1e7 - end_time = ${endTime} - line_search = 'none' - automatic_scaling = true - l_max_its = 20 - l_tol = 8e-3 - nl_forced_its = 1 - nl_max_its = 20 - nl_rel_tol = 5e-04 - nl_abs_tol = 1e-09 -[] - -[Outputs] - print_linear_residuals = false - exodus = false - csv = true - # [fracCSV] - # type = CSV - # sync_times = '100 200 300 400 500 600 700 800 900 - # 1000 2000 3000 4000 5000 6000 7000 8000 9000 - # 1000e1 2000e1 3000e1 4000e1 5000e1 6000e1 7000e1 8000e1 9000e1 - # 1000e2 2000e2 3000e2 4000e2 5000e2 6000e2 7000e2 8000e2 9000e2 - # 1000e3 2000e3 3000e3 4000e3 5000e3 6000e3 7000e3 8000e3 9000e3 - # 1000e4 2000e4 3000e4 4000e4 5000e4 6000e4 7000e4 8000e4 9000e4 - # 1000e5 2000e5 3000e5 4000e5 5000e5 6000e5 7000e5 8000e5 9000e5' - # sync_only = true - # [] -[] diff --git a/examples/gringarten_multiapp/gringartenbrickdomain100by100by10.e b/examples/gringarten_multiapp/gringartenbrickdomain100by100by10.e deleted file mode 100644 index 91ce521e..00000000 Binary files a/examples/gringarten_multiapp/gringartenbrickdomain100by100by10.e and /dev/null differ diff --git a/examples/gringarten_multiapp/mrect1.i b/examples/gringarten_multiapp/mrect1.i index 6d3cee4a..15528d9e 100644 --- a/examples/gringarten_multiapp/mrect1.i +++ b/examples/gringarten_multiapp/mrect1.i @@ -1,11 +1,51 @@ -endTime = 3.16e8 - +endTime = 3e8 +frac_half_space = 20 #50 +n_elem_x = 7 #10 [Mesh] - uniform_refine = 0 - [matrixmesh] - type = FileMeshGenerator # reading mesh from file - file = 'gringartenbrickdomain100by100by10.e' #mesh created in cubit + # [matrixmesh] + # type = FileMeshGenerator # reading mesh from file + # file = 'gringartenbrickdomain100by100by10.e' #mesh created in cubit + # [] + [left_side] + type = GeneratedMeshGenerator + dim = 3 + nx = ${n_elem_x} + xmin = '${fparse -frac_half_space}' + xmax = 0 + bias_x = 0.8 + ny = 50 + ymin = -50 + ymax = 50 + nz = 1 + zmin = -5 + zmax = 5 + [] + [right_side] + type = GeneratedMeshGenerator + dim = 3 + nx = 7 + xmin = 0 + xmax = '${fparse frac_half_space}' + bias_x = 1.2 + ny = 50 + ymin = -50 + ymax = 50 + nz = 1 + zmin = -5 + zmax = 5 + [] + [right_side_rename] + type = RenameBoundaryGenerator + input = right_side + old_boundary = 'top bottom front back' + new_boundary = 'rtop rbottom rfront rback' + [] + [smg] + type = StitchedMeshGenerator + inputs = 'left_side right_side_rename' + clear_stitched_boundary_ids = true + stitch_boundaries_pairs = 'right left' [] [] @@ -39,11 +79,11 @@ endTime = 3.16e8 [Functions] [insitu_pp] type = ParsedFunction - value = '9.81*1000*(3000 - z)' + expression = '9.81*1000*(3000 - z)' [] [insitu_T] type = ParsedFunction - value = '363' + expression = '363' [] [] @@ -84,18 +124,26 @@ endTime = 3.16e8 [] [Preconditioning] - [./superlu] + [superlu] type = SMP full = true petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' petsc_options_value = 'gmres lu superlu_dist' - [../] + [] [] +[Functions] + [dts] + type = ParsedFunction + expression = if(t<1e6,t*t,1e6) + [] +[] [Executioner] type = Transient solve_type = NEWTON - dt = 0.1e7 + # dt = 1e6 + dtmin = 1 + dtmax = 1e6 end_time = ${endTime} line_search = 'none' automatic_scaling = true @@ -105,6 +153,11 @@ endTime = 3.16e8 nl_max_its = 20 nl_rel_tol = 5e-04 nl_abs_tol = 1e-09 + [TimeStepper] + type = FunctionDT + function = dts + min_dt = 2 + [] [] [Outputs] @@ -184,50 +237,43 @@ endTime = 3.16e8 [Transfers] [normal_x_from_fracture] type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app source_variable = normal_dirn_x variable = fracture_normal_x [] [normal_y_from_fracture] type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app source_variable = normal_dirn_y variable = fracture_normal_y [] [normal_z_from_fracture] type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app source_variable = normal_dirn_z variable = fracture_normal_z [] [element_normal_length_to_fracture] type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app + to_multi_app = fracture_app source_variable = element_normal_length variable = enclosing_element_normal_length [] [element_normal_thermal_cond_to_fracture] type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app + to_multi_app = fracture_app source_variable = normal_thermal_conductivity variable = enclosing_element_normal_thermal_cond [] [T_to_fracture] - type = MultiAppInterpolationTransfer - direction = to_multiapp - multi_app = fracture_app + type = MultiAppGeometricInterpolationTransfer + to_multi_app = fracture_app source_variable = matrix_T variable = transferred_matrix_T [] [heat_from_fracture] type = MultiAppReporterTransfer - direction = from_multiapp - multi_app = fracture_app + from_multi_app = fracture_app from_reporters = 'heat_transfer_rate/joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' to_reporters = 'heat_transfer_rate/transferred_joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' [] diff --git a/examples/gringarten_multiapp/mrect1_amr.i b/examples/gringarten_multiapp/mrect1_amr.i deleted file mode 100644 index d0d015b3..00000000 --- a/examples/gringarten_multiapp/mrect1_amr.i +++ /dev/null @@ -1,265 +0,0 @@ -endTime = 3.16e8 - -[Mesh] - uniform_refine = 0 - [generate] - type = GeneratedMeshGenerator - dim = 3 - nx = 10 - xmin = -50 - xmax = 50 - ny = 10 - ymin = -50 - ymax = 50 - nz = 1 - zmin = -5 - zmax = 5 - [] -[] - -[GlobalParams] - PorousFlowDictator = dictator -[] - -[Variables] - [matrix_P] - [] - [matrix_T] - [] -[] - -[ICs] - [matrix_P] - type = FunctionIC - variable = matrix_P - function = insitu_pp - [] - [matrix_T] - type = FunctionIC - variable = matrix_T - function = insitu_T - [] -[] - -[BCs] -[] - -[Functions] - [insitu_pp] - type = ParsedFunction - value = '9.81*1000*(3000 - z)' - [] - [insitu_T] - type = ParsedFunction - value = '363' - [] -[] - -[PorousFlowFullySaturated] - coupling_type = ThermoHydro - porepressure = matrix_P - temperature = matrix_T - fp = water - gravity = '0 0 -9.81' - pressure_unit = Pa -[] - -[FluidProperties] - [water] - type = SimpleFluidProperties - thermal_expansion = 0 - [] -[] - -[Materials] - [porosity] - type = PorousFlowPorosityConst - porosity = 0.1 - [] - [permeability] - type = PorousFlowPermeabilityConst - permeability = '1E-16 0 0 0 1E-16 0 0 0 1E-16' - [] - [internal_energy] - type = PorousFlowMatrixInternalEnergy - density = 2875 - specific_heat_capacity = 825 - [] - [aq_thermal_conductivity] - type = PorousFlowThermalConductivityIdeal - dry_thermal_conductivity = '2.83 0 0 0 2.83 0 0 0 2.83' - [] -[] - -[Preconditioning] - [./superlu] - type = SMP - full = true - petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' - petsc_options_value = 'gmres lu superlu_dist' - [../] -[] - -[Executioner] - type = Transient - solve_type = NEWTON - [TimeStepper] - type = IterationAdaptiveDT - dt = 100 - growth_factor = 2 - optimal_iterations = 20 - [] - dtmin = 1 - dtmax = 0.1e7 - end_time = ${endTime} - line_search = 'none' - # automatic_scaling = true - l_max_its = 20 - l_tol = 8e-3 - nl_forced_its = 1 - nl_max_its = 20 - nl_rel_tol = 5e-04 - nl_abs_tol = 1e-09 -[] - -[Adaptivity] - marker = points - max_h_level = 3 - stop_time = 1000 - [Markers] - [points] - type = ReporterPointMarker - x_coord_name = heat_transfer_rate/x - y_coord_name = heat_transfer_rate/y - z_coord_name = heat_transfer_rate/z - inside = refine - empty = do_nothing - [] - [] -[] - -[Outputs] - print_linear_residuals = false -[] - -[DiracKernels] - [heat_from_fracture] - type = ReporterPointSource - variable = matrix_T - value_name = heat_transfer_rate/transferred_joules_per_s - x_coord_name = heat_transfer_rate/x - y_coord_name = heat_transfer_rate/y - z_coord_name = heat_transfer_rate/z - [] -[] - -[Reporters] - [heat_transfer_rate] - type = ConstantReporter - real_vector_names = 'transferred_joules_per_s x y z' - real_vector_values = '10000; 0; 0; 0' - outputs = none - [] -[] - -[AuxVariables] - [normal_thermal_conductivity] - family = MONOMIAL - order = CONSTANT - [] - [fracture_normal_x] - family = MONOMIAL - order = CONSTANT - initial_condition = 0 - [] - [fracture_normal_y] - family = MONOMIAL - order = CONSTANT - initial_condition = 1 - [] - [fracture_normal_z] - family = MONOMIAL - order = CONSTANT - initial_condition = 0 - [] - [element_normal_length] - family = MONOMIAL - order = CONSTANT - [] - [fracDensity_AMR] - [] -[] - -[AuxKernels] - [normal_thermal_conductivity_auxk] - type = ConstantAux - variable = normal_thermal_conductivity - value = 5 - [] - [element_normal_length_auxk] - type = PorousFlowElementLength - variable = element_normal_length - direction = 'fracture_normal_x fracture_normal_y fracture_normal_z' - [] -[] - -[MultiApps] - [fracture_app] - type = TransientMultiApp - input_files = frect1.i - execute_on = TIMESTEP_BEGIN - sub_cycling = true - [] -[] - -[Transfers] - [normal_x_from_fracture] - type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app - source_variable = normal_dirn_x - variable = fracture_normal_x - [] - [normal_y_from_fracture] - type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app - source_variable = normal_dirn_y - variable = fracture_normal_y - [] - [normal_z_from_fracture] - type = MultiAppNearestNodeTransfer - direction = from_multiapp - multi_app = fracture_app - source_variable = normal_dirn_z - variable = fracture_normal_z - [] - [element_normal_length_to_fracture] - type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app - source_variable = element_normal_length - variable = enclosing_element_normal_length - [] - [element_normal_thermal_cond_to_fracture] - type = MultiAppNearestNodeTransfer - direction = to_multiapp - multi_app = fracture_app - source_variable = normal_thermal_conductivity - variable = enclosing_element_normal_thermal_cond - [] - [T_to_fracture] - type = MultiAppInterpolationTransfer - direction = to_multiapp - multi_app = fracture_app - source_variable = matrix_T - variable = transferred_matrix_T - [] - [heat_from_fracture] - type = MultiAppReporterTransfer - direction = from_multiapp - multi_app = fracture_app - from_reporters = 'heat_transfer_rate/joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' - to_reporters = 'heat_transfer_rate/transferred_joules_per_s heat_transfer_rate/x heat_transfer_rate/y heat_transfer_rate/z' - [] -[] diff --git a/examples/gringarten_multiapp/plot_compare.py b/examples/gringarten_multiapp/plot_compare.py index 72f9c8ed..923e5fcc 100755 --- a/examples/gringarten_multiapp/plot_compare.py +++ b/examples/gringarten_multiapp/plot_compare.py @@ -1,31 +1,34 @@ #!/usr/bin/env python3 -#* This file is part of the MOOSE framework -#* https://www.mooseframework.org -#* -#* All rights reserved, see COPYRIGHT for full restrictions -#* https://github.com/idaholab/moose/blob/master/COPYRIGHT -#* -#* Licensed under LGPL 2.1, please see LICENSE for details -#* https://www.gnu.org/licenses/lgpl-2.1.html +# * This file is part of the MOOSE framework +# * https://www.mooseframework.org +# * +# * All rights reserved, see COPYRIGHT for full restrictions +# * https://github.com/idaholab/moose/blob/master/COPYRIGHT +# * +# * Licensed under LGPL 2.1, please see LICENSE for details +# * https://www.gnu.org/licenses/lgpl-2.1.html import pandas as pd -import os -import sys import numpy as np import matplotlib.pyplot as plt import math -#------------------------------------------------------------------------------- -df = pd.read_csv('mrect1_out_fracture_app0.csv') -print('fracture dataframe description: ',df.columns.values) -['time' 'J_out' 'P_in' 'P_out' 'TK_in' 'TK_out' 'kg_out' 'kg_per_s'] -df_amr = pd.read_csv('mrect1_amr_out_fracture_app0.csv') -#analytic solution from Koenradd and Gringarten paper, 1975 -timevector = np.logspace(1, 9, 50, base=10) -print(timevector) - -#analytic Gringarten solution 1975 from Beckers, Koenraad -Trock = 363-273.15 -Tinj = 303-273.15 +from mpmath import * + + +# ------------------------------------------------------------------------------- +df = pd.read_csv("mrect1_out_fracture_app0.csv") +print("fracture dataframe description: ", df.columns.values) +["time" "J_out" "P_in" "P_out" "TK_in" "TK_out" "kg_out" "kg_per_s"] +df_50 = pd.read_csv("mrect1_out_fracture_app0_d50.csv") +# analytic solution from Koenradd and Gringarten paper, 1975 +timevector = np.logspace(6, 11, 50, base=10) + +# analytic solution for single fracture Gringarten 1975 eqn A18 +# from Beckers, Koenraad + +# terms for infinite and finite space fractures +Trock = 363 - 273.15 +Tinj = 303 - 273.15 rhowater = 1000 cpwater = 4150 fracnumb = 1 @@ -35,21 +38,86 @@ fracheight = 10 rhorock = 2875 cprock = 825 -mtot = 25/25/10 -Q = mtot/rhowater/fracnumb/fracwidth -td = (rhowater*cpwater)*(rhowater*cpwater)/(4*krock*rhorock*cprock)*(Q/fracheight)*(Q/fracheight)*timevector -Twater=[] -for i in range(len(timevector)): - temp=1-math.erf(1/(2*math.sqrt(td[i]))) - Twater.append(Trock-temp*(Trock-Tinj)) - - -#------------------------------------------------------------------------------- -fig1=plt.figure() -plt.plot(timevector/3600/24/365,Twater,'b', linestyle = '-',marker = 'o',fillstyle='none',label='Gringarten 1975') -plt.plot(df['time']/3600/24/365,df['TK_out']-273.15,'r', linestyle = '--',label='Falcon') -plt.plot(df_amr['time']/3600/24/365,df_amr['TK_out']-273.15,'g', linestyle = '-.',label='Falcon amr') -plt.ylim([50, 100]) +mtot = 25 / 25 / 10 +Q = mtot / rhowater / fracnumb / fracwidth +# terms needed for finite spaced fractures. Gringarten 1975 eqn A19 +xe = 20 # 1/2 spacing between fractures +z = fracheight + + +# nonDimensionalized Time eqn 10 +td = ( + (rhowater * cpwater) + * (rhowater * cpwater) + / (4 * krock * rhorock * cprock) # where does the 4 come from? Its not in eqn 10 + * (Q / fracheight) + * (Q / fracheight) + * timevector +) +# Solution to eqn A18 -- infinite space between fractures +Twater = [] +Twd_inf = [] +for i in range(len(td)): + # eqn A18 solving for nonDimensionalized temperature + temp = 1 - math.erf(1 / (2 * math.sqrt(td[i]))) + Twd_inf.append(temp) + # convert to Temperature + Twater.append(Trock - temp * (Trock - Tinj)) + +# Solution to eqn A19 -- finite space between fractures set above by xe +# simplified inverse laplace transform +mp.dps = 15 +top = rhowater * cpwater * Q * xe +bottom = 2 * krock * fracheight +zd = z / fracheight +Twd_tilde = lambda s: 1 / s * exp(-zd * sqrt(s) * tanh(top / bottom * sqrt(s))) +Twater_finite = [] +for i in range(len(td)): + temp = invertlaplace(Twd_tilde, td[i], method="talbot") + # convert to Temperature + Twater_finite.append(Trock - temp * (Trock - Tinj)) + +# ------------------------------------------------------------------------------- +fig1 = plt.figure() +plt.plot( + timevector / 3600 / 24 / 365, + Twater, + "b", + linestyle="-", + linewidth=1, + marker="o", + alpha=0.5, + fillstyle="none", + label="Gringarten xe=inf", +) +plt.plot( + timevector / 3600 / 24 / 365, + Twater_finite, + "k", + linestyle="-", + linewidth=1, + marker="^", + alpha=0.5, + fillstyle="none", + label="Gringarten xe=20m", +) +plt.plot( + df_50["time"] / 3600 / 24 / 365, + df_50["TK_out"] - 273.15, + "r", + linestyle="--", + linewidth=2, + label="Falcon xe=50m", +) +plt.plot( + df["time"] / 3600 / 24 / 365, + df["TK_out"] - 273.15, + "g", + linestyle="-.", + linewidth=2, + label="Falcon xe=20m", +) +plt.ylim([30, 100]) plt.xlim([0, 10]) plt.ylabel("T (degC)") plt.xlabel("Time (year)") diff --git a/examples/gringarten_multiapp/production_single_fracture_amr.xyz b/examples/gringarten_multiapp/production_single_fracture_amr.xyz deleted file mode 100644 index 83008e41..00000000 --- a/examples/gringarten_multiapp/production_single_fracture_amr.xyz +++ /dev/null @@ -1 +0,0 @@ -0.01 0.1 50 0 diff --git a/examples/gringarten_multiapp/rect_100_10_medium.e b/examples/gringarten_multiapp/rect_100_10_medium.e deleted file mode 100644 index 2a3ae8c7..00000000 Binary files a/examples/gringarten_multiapp/rect_100_10_medium.e and /dev/null differ diff --git a/examples/gringarten_multiapp/tests b/examples/gringarten_multiapp/tests index da863bd3..3a6aa910 100644 --- a/examples/gringarten_multiapp/tests +++ b/examples/gringarten_multiapp/tests @@ -1,11 +1,10 @@ [Tests] [shortTest] - method = '!dbg' type = CSVDiff input = 'mrect1.i' csvdiff = 'mrect1_out_fracture_app0.csv' - cli_args = 'endTime=6e7' - requirement = 'Run Gringarten multiapp for a few timesteps.' - skip = 'long time' + cli_args = 'endTime=6e6' +requirement = 'The system shall match the Gringarten analytic solution for a set of parallel fracturs seperated by a finite distance using a multi-app to transfer heat between a 2D fracture and 3D matrix.' + heavy = true [../] []