Skip to content
Snippets Groups Projects
Commit 973966c5 authored by David's avatar David
Browse files

overhaul TPR layered soil

parent 37163fe9
No related branches found
No related tags found
No related merge requests found
#!/usr/bin/python3
"""Layered soil simulation.
""" TP-R Layered soil simulation.
This program sets up an LDD simulation
"""
import dolfin as df
# import mshr
# import numpy as np
import sympy as sym
# import typing as tp
import functools as ft
# import domainPatch as dp
import LDDsimulation as ldd
import helpers as hlp
import datetime
import os
import pandas as pd
# check if output directory exists
# init sympy session
sym.init_printing()
# PREREQUISITS ###############################################################
# check if output directory "./output" exists. This will be used in
# the generation of the output string.
if not os.path.exists('./output'):
os.mkdir('./output')
print("Directory ", './output', " created ")
......@@ -25,42 +26,51 @@ else:
print("Directory ", './output', " already exists. Will use as output \
directory")
date = datetime.datetime.now()
datestr = date.strftime("%Y-%m-%d")
# init sympy session
sym.init_printing()
# solver_tol = 6E-7
use_case = "TP-R-layered-soil-all-params-set-one-new-lambda"
# Name of the usecase that will be printed during simulation.
use_case = "TP-R-layered-soil-all-params-one"
# The name of this very file. Needed for creating log output.
thisfile = "TP-R-layered_soil-all-params-one.py"
max_iter_num = 500
# GENERAL SOLVER CONFIG ######################################################
# maximal iteration per timestep
max_iter_num = 300
FEM_Lagrange_degree = 1
# GRID AND MESH STUDY SPECIFICATIONS #########################################
mesh_study = False
resolutions = {
# 1: 1e-7, # h=2
# 2: 2e-5, # h=1.1180
# 4: 1e-6, # h=0.5590
# 8: 1e-6, # h=0.2814
# 16: 1e-6, # h=0.1412
32: 1e-6,
# 64: 5e-7,
# 128: 5e-7
# 1: 1e-6,
# 2: 1e-6,
# 4: 1e-6,
# 8: 1e-6,
16: 5e-6,
# 32: 5e-6,
# 64: 2e-6,
# 128: 1e-6,
# 256: 1e-6,
}
# GRID #######################
# mesh_resolution = 20
timestep_size = 0.001
number_of_timesteps = 5
plot_timestep_every = 4
# decide how many timesteps you want analysed. Analysed means, that we write
# out subsequent errors of the L-iteration within the timestep.
number_of_timesteps_to_analyse = 6
# starttimes gives a list of starttimes to run the simulation from.
# The list is looped over and a simulation is run with t_0 as initial time
# for each element t_0 in starttimes.
starttimes = [0.0]
timestep_size = 0.001
number_of_timesteps = 20
Lw = 0.025 # /timestep_size
Lnw = Lw
# LDD scheme parameters ######################################################
Lw1 = 0.025 # /timestep_size
Lnw1 = Lw1
Lw2 = 0.025 # /timestep_size
Lnw2 = Lw2
Lw3 = 0.025 # /timestep_size
Lnw3 = Lw3
Lw4 = 0.025 # /timestep_size
Lnw4 = Lw4
lambda12_w = 40
lambda12_nw = 40
......@@ -69,31 +79,42 @@ lambda23_nw = 40
lambda34_w = 40
lambda34_nw = 40
include_gravity = True
include_gravity = False
debugflag = False
analyse_condition = True
if mesh_study:
output_string = "./output/{}-{}_timesteps{}_P{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree
)
else:
for tol in resolutions.values():
solver_tol = tol
output_string = "./output/{}-{}_timesteps{}_P{}_solver_tol{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree, solver_tol
)
# I/O CONFIG #################################################################
# when number_of_timesteps is high, it might take a long time to write all
# timesteps to disk. Therefore, you can choose to only write data of every
# plot_timestep_every timestep to disk.
plot_timestep_every = 4
# Decide how many timesteps you want analysed. Analysed means, that
# subsequent errors of the L-iteration within the timestep are written out.
number_of_timesteps_to_analyse = 5
# toggle what should be written to files
# fine grained control over data to be written to disk in the mesh study case
# as well as for a regular simuation for a fixed grid.
if mesh_study:
write_to_file = {
# output the relative errornorm (integration in space) w.r.t. an exact
# solution for each timestep into a csv file.
'space_errornorms': True,
# save the mesh and marker functions to disk
'meshes_and_markers': True,
# save xdmf/h5 data for each LDD iteration for timesteps determined by
# number_of_timesteps_to_analyse. I/O intensive!
'L_iterations_per_timestep': False,
# save solution to xdmf/h5.
'solutions': True,
# save absolute differences w.r.t an exact solution to xdmf/h5 file
# to monitor where on the domains errors happen
'absolute_differences': True,
# analyise condition numbers for timesteps determined by
# number_of_timesteps_to_analyse and save them over time to csv.
'condition_numbers': analyse_condition,
# output subsequent iteration errors measured in L^2 to csv for
# timesteps determined by number_of_timesteps_to_analyse.
# Usefull to monitor convergence of the acutal LDD solver.
'subsequent_errors': True
}
else:
......@@ -107,7 +128,20 @@ else:
'subsequent_errors': True
}
# OUTPUT FILE STRING #########################################################
if mesh_study:
output_string = "./output/{}-{}_timesteps{}_P{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree
)
else:
for tol in resolutions.values():
solver_tol = tol
output_string = "./output/{}-{}_timesteps{}_P{}_solver_tol{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree, solver_tol
)
# DOMAIN AND INTERFACE #######################################################
# global domain
subdomain0_vertices = [
df.Point(-1.0, -1.0),
......@@ -257,6 +291,8 @@ outer_boundary_def_points = {
4: subdomain4_outer_boundary_verts
}
# MODEL CONFIGURATION #########################################################
isRichards = {
1: True,
2: True,
......@@ -273,53 +309,51 @@ isRichards = {
# Dict of the form: { subdom_num : viscosity }
viscosity = {
1: {'wetting': 1,
'nonwetting': 1},
2: {'wetting': 1,
'nonwetting': 1},
3: {'wetting': 1,
'nonwetting': 1},
4: {'wetting': 1,
'nonwetting': 1},
1: {'wetting': 1.0,
'nonwetting': 1.0},
2: {'wetting': 1.0,
'nonwetting': 1.0},
3: {'wetting': 1.0,
'nonwetting': 1.0},
4: {'wetting': 1.0,
'nonwetting': 1.0},
}
# Dict of the form: { subdom_num : density }
densities = {
1: {'wetting': 1, # 997
'nonwetting': 1}, # 1.225}},
2: {'wetting': 1, # 997
'nonwetting': 1}, # 1.225}},
3: {'wetting': 1, # 997
'nonwetting': 1}, # 1.225}},
4: {'wetting': 1, # 997
'nonwetting': 1}, # 1.225}}
1: {'wetting': 1.0, # 997
'nonwetting': 1.0}, # 1.225}},
2: {'wetting': 1.0, # 997
'nonwetting': 1.0}, # 1.225}},
3: {'wetting': 1.0, # 997
'nonwetting': 1.0}, # 1.225}},
4: {'wetting': 1.0, # 997
'nonwetting': 1.0}, # 1.225}}
}
gravity_acceleration = 1
gravity_acceleration = 1.0
# porosities taken from
# https://www.geotechdata.info/parameter/soil-porosity.html
# Dict of the form: { subdom_num : porosity }
porosity = {
1: 1, # 0.2, # Clayey gravels, clayey sandy gravels
2: 1, # 0.22, # Silty gravels, silty sandy gravels
3: 1, # 0.37, # Clayey sands
4: 1, # 0.2 # Silty or sandy clay
1: 1.0, #0.37, # 0.2, # Clayey gravels, clayey sandy gravels
2: 1.0, #0.22, # 0.22, # Silty gravels, silty sandy gravels
3: 1.0, #0.2, # 0.37, # Clayey sands
4: 1.0, #0.22, # 0.2 # Silty or sandy clay
}
# subdom_num : subdomain L for L-scheme
L = {
1: {'wetting': Lw,
'nonwetting': Lnw},
2: {'wetting': Lw,
'nonwetting': Lnw},
3: {'wetting': Lw,
'nonwetting': Lnw},
4: {'wetting': Lw,
'nonwetting': Lnw}
1: {'wetting': Lw1,
'nonwetting': Lnw1},
2: {'wetting': Lw2,
'nonwetting': Lnw2},
3: {'wetting': Lw3,
'nonwetting': Lnw3},
4: {'wetting': Lw4,
'nonwetting': Lnw4}
}
# subdom_num : lambda parameter for the L-scheme
# interface_num : lambda parameter for the L-scheme on that interface.
# Note that interfaces are numbered starting from 0, because
# adjacent_subdomains is a list and not a dict. Historic fuckup, I know
......@@ -332,11 +366,13 @@ lambda_param = {
'nonwetting': lambda34_nw},
}
# after Lewis, see pdf file
intrinsic_permeability = {
1: 1,
2: 1,
3: 1,
4: 1
1: 1.0, #0.01, # sand
2: 1.0, #0.01, # sand, there is a range
3: 1.0, #0.01, #10e-2, # clay has a range
4: 1.0, #0.01, #10e-3
}
......@@ -354,19 +390,46 @@ def rel_perm1nw(s):
# relative permeabilty functions on subdomain 2
def rel_perm2w(s):
# relative permeabilty wetting on subdomain2
return intrinsic_permeability[2]*s**3
return intrinsic_permeability[2]*s**2
def rel_perm2nw(s):
# relative permeabilty nonwetting on subdomain2
return intrinsic_permeability[2]*(1-s)**3
return intrinsic_permeability[2]*(1-s)**2
# relative permeabilty functions on subdomain 3
def rel_perm3w(s):
# relative permeabilty wetting on subdomain3
return intrinsic_permeability[3]*s**3
def rel_perm3nw(s):
# relative permeabilty nonwetting on subdomain3
return intrinsic_permeability[3]*(1-s)**3
# relative permeabilty functions on subdomain 4
def rel_perm4w(s):
# relative permeabilty wetting on subdomain4
return intrinsic_permeability[4]*s**3
def rel_perm4nw(s):
# relative permeabilty nonwetting on subdomain4
return intrinsic_permeability[4]*(1-s)**3
_rel_perm1w = ft.partial(rel_perm1w)
_rel_perm1nw = ft.partial(rel_perm1nw)
_rel_perm2w = ft.partial(rel_perm2w)
_rel_perm2nw = ft.partial(rel_perm2nw)
_rel_perm3w = ft.partial(rel_perm3w)
_rel_perm3nw = ft.partial(rel_perm3nw)
_rel_perm4w = ft.partial(rel_perm4w)
_rel_perm4nw = ft.partial(rel_perm4nw)
subdomain1_rel_perm = {
'wetting': _rel_perm1w,
'nonwetting': _rel_perm1nw
......@@ -377,18 +440,28 @@ subdomain2_rel_perm = {
'nonwetting': _rel_perm2nw
}
# _rel_perm3 = ft.partial(rel_perm2)
# subdomain3_rel_perm = subdomain2_rel_perm.copy()
#
# _rel_perm4 = ft.partial(rel_perm1)
# subdomain4_rel_perm = subdomain1_rel_perm.copy()
subdomain3_rel_perm = {
'wetting': _rel_perm3w,
'nonwetting': _rel_perm3nw
}
subdomain4_rel_perm = {
'wetting': _rel_perm4w,
'nonwetting': _rel_perm4nw
}
# dictionary of relative permeabilties on all domains.
# relative_permeability = {
# 1: subdomain1_rel_perm,
# 2: subdomain1_rel_perm,
# 3: subdomain2_rel_perm,
# 4: subdomain2_rel_perm
# }
relative_permeability = {
1: subdomain1_rel_perm,
2: subdomain1_rel_perm,
3: subdomain2_rel_perm,
4: subdomain2_rel_perm
2: subdomain2_rel_perm,
3: subdomain3_rel_perm,
4: subdomain4_rel_perm
}
......@@ -404,22 +477,48 @@ def rel_perm1nw_prime(s):
return -1*intrinsic_permeability[1]*2*(1-s)
# definition of the derivatives of the relative permeabilities
# relative permeabilty functions on subdomain 1
def rel_perm2w_prime(s):
# relative permeabilty on subdomain1
return intrinsic_permeability[2]*3*s**2
# relative permeabilty on subdomain2
return intrinsic_permeability[2]*2*s
def rel_perm2nw_prime(s):
# relative permeabilty on subdomain1
return -1*intrinsic_permeability[2]*3*(1-s)**2
# relative permeabilty on subdomain2
return -1*intrinsic_permeability[2]*2*(1-s)
# definition of the derivatives of the relative permeabilities
# relative permeabilty functions on subdomain 3
def rel_perm3w_prime(s):
# relative permeabilty on subdomain3
return intrinsic_permeability[3]*3*s**2
def rel_perm3nw_prime(s):
# relative permeabilty on subdomain3
return -1*intrinsic_permeability[3]*3*(1-s)**2
# definition of the derivatives of the relative permeabilities
# relative permeabilty functions on subdomain 4
def rel_perm4w_prime(s):
# relative permeabilty on subdomain4
return intrinsic_permeability[4]*3*s**2
def rel_perm4nw_prime(s):
# relative permeabilty on subdomain4
return -1*intrinsic_permeability[4]*3*(1-s)**2
_rel_perm1w_prime = ft.partial(rel_perm1w_prime)
_rel_perm1nw_prime = ft.partial(rel_perm1nw_prime)
_rel_perm2w_prime = ft.partial(rel_perm2w_prime)
_rel_perm2nw_prime = ft.partial(rel_perm2nw_prime)
_rel_perm3w_prime = ft.partial(rel_perm3w_prime)
_rel_perm3nw_prime = ft.partial(rel_perm3nw_prime)
_rel_perm4w_prime = ft.partial(rel_perm4w_prime)
_rel_perm4nw_prime = ft.partial(rel_perm4nw_prime)
subdomain1_rel_perm_prime = {
'wetting': _rel_perm1w_prime,
......@@ -432,15 +531,32 @@ subdomain2_rel_perm_prime = {
'nonwetting': _rel_perm2nw_prime
}
subdomain3_rel_perm_prime = {
'wetting': _rel_perm3w_prime,
'nonwetting': _rel_perm3nw_prime
}
subdomain4_rel_perm_prime = {
'wetting': _rel_perm4w_prime,
'nonwetting': _rel_perm4nw_prime
}
# dictionary of relative permeabilties on all domains.
# ka_prime = {
# 1: subdomain1_rel_perm_prime,
# 2: subdomain1_rel_perm_prime,
# 3: subdomain2_rel_perm_prime,
# 4: subdomain2_rel_perm_prime
# }
ka_prime = {
1: subdomain1_rel_perm_prime,
2: subdomain1_rel_perm_prime,
3: subdomain2_rel_perm_prime,
4: subdomain2_rel_perm_prime
2: subdomain2_rel_perm_prime,
3: subdomain3_rel_perm_prime,
4: subdomain4_rel_perm_prime
}
# S-pc-relation ship. We use the van Genuchten approach, i.e.
# pc = 1/alpha*(S^{-1/m} -1)^1/n, where we set alpha = 0, assume m = 1-1/n
# (see Helmig) and assume that residual saturation is Sw
......@@ -500,17 +616,18 @@ sat_pressure_relationship = {
}
#############################################
###############################################################################
# Manufacture source expressions with sympy #
#############################################
###############################################################################
x, y = sym.symbols('x[0], x[1]') # needed by UFL
t = sym.symbols('t', positive=True)
p_e_sym_2patch = {
1: {'wetting': -7 - (1+t*t)*(1 + x*x + y*y),
1: {'wetting': (-7.0 - (1.0 + t*t)*(1.0 + x*x + y*y)),
'nonwetting': 0.0*t}, # -1-t*(1.1 + y + x**2)**2},
2: {'wetting': -7.0 - (1.0 + t*t)*(1.0 + x*x),
'nonwetting': (-1-t*(1.1 + x**2)**2 - sym.sqrt(5+t**2))*y**2},
'nonwetting': (-2-t*(1.1+y + x**2))*y**2},
# 'nonwetting': (-1-t*(1.1 + x**2)**2 - sym.sqrt(5+t**2))*y**2},
}
p_e_sym = {
......@@ -583,6 +700,7 @@ source_expression = exact_solution_example['source']
exact_solution = exact_solution_example['exact_solution']
initial_condition = exact_solution_example['initial_condition']
# BOUNDARY CONDITIONS #########################################################
# Dictionary of dirichlet boundary conditions.
dirichletBC = dict()
# similarly to the outer boundary dictionary, if a patch has no outer boundary
......@@ -612,6 +730,7 @@ for subdomain in isRichards.keys():
{outer_boundary_ind: exact_solution[subdomain]}
)
# LOG FILE OUTPUT #############################################################
# read this file and print it to std out. This way the simulation can produce a
# log file with ./TP-R-layered_soil.py | tee simulation.log
......
#!/usr/bin/python3
"""Layered soil simulation.
""" TP-R Layered soil simulation.
This program sets up an LDD simulation
"""
import dolfin as df
# import mshr
# import numpy as np
import sympy as sym
# import typing as tp
import functools as ft
# import domainPatch as dp
import LDDsimulation as ldd
import helpers as hlp
import datetime
import os
import pandas as pd
# check if output directory exists
# init sympy session
sym.init_printing()
# PREREQUISITS ###############################################################
# check if output directory "./output" exists. This will be used in
# the generation of the output string.
if not os.path.exists('./output'):
os.mkdir('./output')
print("Directory ", './output', " created ")
......@@ -25,40 +26,43 @@ else:
print("Directory ", './output', " already exists. Will use as output \
directory")
date = datetime.datetime.now()
datestr = date.strftime("%Y-%m-%d")
# init sympy session
sym.init_printing()
# solver_tol = 6E-7
# Name of the usecase that will be printed during simulation.
use_case = "TP-R-layered-soil-realistic"
# name of this very file. Needed for log output.
# The name of this very file. Needed for creating log output.
thisfile = "TP-R-layered_soil.py"
# GENERAL SOLVER CONFIG ######################################################
# maximal iteration per timestep
max_iter_num = 300
FEM_Lagrange_degree = 1
# GRID AND MESH STUDY SPECIFICATIONS #########################################
mesh_study = False
resolutions = {
# 1: 2e-6, # h=2
# 2: 2e-6, # h=1.1180
# 4: 2e-6, # h=0.5590
# 8: 2e-6, # h=0.2814
# 16: 2e-6, # h=0.1412
32: 2e-6,
# 64: 2e-6,
# 128: 2e-6
# 1: 1e-6,
# 2: 1e-6,
# 4: 1e-6,
# 8: 1e-6,
# 16: 5e-6,
# 32: 5e-6,
64: 2e-6,
# 128: 1e-6,
# 256: 1e-6,
}
# GRID #######################
# mesh_resolution = 20
timestep_size = 0.001
number_of_timesteps = 1000
plot_timestep_every = 4
# decide how many timesteps you want analysed. Analysed means, that we write
# out subsequent errors of the L-iteration within the timestep.
number_of_timesteps_to_analyse = 4
# starttimes gives a list of starttimes to run the simulation from.
# The list is looped over and a simulation is run with t_0 as initial time
# for each element t_0 in starttimes.
starttimes = [0.0]
timestep_size = 0.001
number_of_timesteps = 20
# LDD scheme parameters ######################################################
Lw1 = 0.025 # /timestep_size
Lnw1 = Lw1
Lw2 = 0.025 # /timestep_size
......@@ -79,27 +83,38 @@ include_gravity = False
debugflag = False
analyse_condition = True
if mesh_study:
output_string = "./output/{}-{}_timesteps{}_P{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree
)
else:
for tol in resolutions.values():
solver_tol = tol
output_string = "./output/{}-{}_timesteps{}_P{}_solver_tol{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree, solver_tol
)
# I/O CONFIG #################################################################
# when number_of_timesteps is high, it might take a long time to write all
# timesteps to disk. Therefore, you can choose to only write data of every
# plot_timestep_every timestep to disk.
plot_timestep_every = 4
# Decide how many timesteps you want analysed. Analysed means, that
# subsequent errors of the L-iteration within the timestep are written out.
number_of_timesteps_to_analyse = 5
# toggle what should be written to files
# fine grained control over data to be written to disk in the mesh study case
# as well as for a regular simuation for a fixed grid.
if mesh_study:
write_to_file = {
# output the relative errornorm (integration in space) w.r.t. an exact
# solution for each timestep into a csv file.
'space_errornorms': True,
# save the mesh and marker functions to disk
'meshes_and_markers': True,
# save xdmf/h5 data for each LDD iteration for timesteps determined by
# number_of_timesteps_to_analyse. I/O intensive!
'L_iterations_per_timestep': False,
# save solution to xdmf/h5.
'solutions': True,
# save absolute differences w.r.t an exact solution to xdmf/h5 file
# to monitor where on the domains errors happen
'absolute_differences': True,
# analyise condition numbers for timesteps determined by
# number_of_timesteps_to_analyse and save them over time to csv.
'condition_numbers': analyse_condition,
# output subsequent iteration errors measured in L^2 to csv for
# timesteps determined by number_of_timesteps_to_analyse.
# Usefull to monitor convergence of the acutal LDD solver.
'subsequent_errors': True
}
else:
......@@ -113,7 +128,20 @@ else:
'subsequent_errors': True
}
# OUTPUT FILE STRING #########################################################
if mesh_study:
output_string = "./output/{}-{}_timesteps{}_P{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree
)
else:
for tol in resolutions.values():
solver_tol = tol
output_string = "./output/{}-{}_timesteps{}_P{}_solver_tol{}".format(
datestr, use_case, number_of_timesteps, FEM_Lagrange_degree, solver_tol
)
# DOMAIN AND INTERFACE #######################################################
# global domain
subdomain0_vertices = [
df.Point(-1.0, -1.0),
......@@ -263,6 +291,8 @@ outer_boundary_def_points = {
4: subdomain4_outer_boundary_verts
}
# MODEL CONFIGURATION #########################################################
isRichards = {
1: True,
2: True,
......@@ -346,23 +376,14 @@ lambda_param = {
2: {'wetting': lambda34_w,
'nonwetting': lambda34_nw},
}
# lambda_param = {
# 1: {'wetting': lambda_w,
# 'nonwetting': lambda_nw},
# 2: {'wetting': lambda_w,
# 'nonwetting': lambda_nw},
# 3: {'wetting': lambda_w,
# 'nonwetting': lambda_nw},
# 4: {'wetting': lambda_w,
# 'nonwetting': lambda_nw},
# }
# after Lewis, see pdf file
intrinsic_permeability = {
1: 0.1, # sand
2: 0.1, # sand, there is a range
3: 0.001, #10e-2, # clay has a range
4: 0.001, #10e-3
1: 0.01, # sand
2: 0.01, # sand, there is a range
3: 0.01, #10e-2, # clay has a range
4: 0.01, #10e-3
}
......@@ -606,17 +627,17 @@ sat_pressure_relationship = {
}
#############################################
###############################################################################
# Manufacture source expressions with sympy #
#############################################
###############################################################################
x, y = sym.symbols('x[0], x[1]') # needed by UFL
t = sym.symbols('t', positive=True)
p_e_sym_2patch = {
1: {'wetting': -6 - (1+t*t)*(1 + x*x + y*y),
1: {'wetting': (-7.0 - (1.0 + t*t)*(1.0 + x*x + y*y)),
'nonwetting': 0.0*t}, # -1-t*(1.1 + y + x**2)**2},
2: {'wetting': -6.0 - (1.0 + t*t)*(1.0 + x*x),
'nonwetting': (-1-t*(1.1+y + x**2))*y**2},
2: {'wetting': -7.0 - (1.0 + t*t)*(1.0 + x*x),
'nonwetting': (-2-t*(1.1+y + x**2))*y**2},
# 'nonwetting': (-1-t*(1.1 + x**2)**2 - sym.sqrt(5+t**2))*y**2},
}
......@@ -690,6 +711,7 @@ source_expression = exact_solution_example['source']
exact_solution = exact_solution_example['exact_solution']
initial_condition = exact_solution_example['initial_condition']
# BOUNDARY CONDITIONS #########################################################
# Dictionary of dirichlet boundary conditions.
dirichletBC = dict()
# similarly to the outer boundary dictionary, if a patch has no outer boundary
......
......@@ -464,6 +464,7 @@ source_expression = exact_solution_example['source']
exact_solution = exact_solution_example['exact_solution']
initial_condition = exact_solution_example['initial_condition']
# BOUNDARY CONDITIONS #########################################################
# Dictionary of dirichlet boundary conditions.
dirichletBC = dict()
# similarly to the outer boundary dictionary, if a patch has no outer boundary
......@@ -478,7 +479,6 @@ dirichletBC = dict()
# return the actual expression needed for the dirichlet condition for both
# phases if present.
# BOUNDARY CONDITIONS #########################################################
# subdomain index: {outer boudary part index: {phase: expression}}
for subdomain in isRichards.keys():
# subdomain can have no outer boundary
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment