diff --git a/src/input_DOUAR_v0.2.txt b/src/input_DOUAR_v0.2.txt new file mode 100644 index 0000000000000000000000000000000000000000..d5046a0ecf77301d8f0070b2aa21d0d363636224 --- /dev/null +++ b/src/input_DOUAR_v0.2.txt @@ -0,0 +1,805 @@ +!------------------------------------------------------------------------------| +!------------------------------------------------------------------------------| +! | +! 8888888b. .d88888b. 888 888 d8888 8888888b. | +! 888 "Y88b d88P" "Y88b 888 888 d88888 888 Y88b | +! 888 888 888 888 888 888 d88P888 888 888 | +! 888 888 888 888 888 888 d88P 888 888 d88P | +! 888 888 888 888 888 888 d88P 888 8888888P" | +! 888 888 888 888 888 888 d88P 888 888 T88b | +! 888 .d88P Y88b. .d88P Y88b. .d88P d8888888888 888 T88b | +! 8888888P" "Y88888P" "Y88888P" d88P 888 888 T88b | +! | +! Input file for DOUAR-WSMP version 0.2 - 23/03/2011 | +!------------------------------------------------------------------------------| +!------------------------------------------------------------------------------| + + +This file is read by two subroutines in the code: +- read_controlling_parameters +- read_input_file +The indications between [] indicate the type of the read parameter. It can be +an integer [int], a real*8 [dp], a character chain [char] or a boolean [bool]. + + +CONTROLLING PARAMETERS +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]debug is a parameter that switches on/off various prints and outputs +(the level of printing for error, warning and messages of the solver id%cntl(4) +is set to the debug value). +if debug is zero, no debugging +if debug is equal to one, this triggers the terminal display of some key +parameters +if debug is equal to two, same as debug equal to one, and the code produces +various vtk files in the DEBUG subdirectories (surfaces, olsf, ...), as well as +cross sections for all nonlinear iterations. Careful, memory consuming! + + debug = 1 + +[bool]doDoRuRe is a flag that triggers the production of output files needed to +produce the DoRuRes. DoRuRe stands for 'Douar Run Report'. + + doDoRuRe = F + +[bool]compute_qp_gram triggers the production of qpgrams for every grid. + + compute_qpgram = F + +[bool]compute_reaction_forces toggles on/off the reaction forces computations. + + compute_reaction_forces = F + + +RESTART +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]irestart is a restart flag; if irestart is not 0, the run will restart from +an output file given by [char]restartfile and at step. + + irestart = 0 + + restartfile = OUT/time_0000.bin + + +TIMESTEPPING +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[dp]dt is the time step length (if dt is negative, courant conidition is used +and automatic time stepping is turned on) + + dt = 1.25d-3 + +[int]nstep is the number of time steps + + nstep = 50 + +[dp]courant is only used when dt is negative; it determines the size of the time +step from the maximum value of the velocity field amplitude. The time step s +the product of courant by the ratio of the smallest leaf size by the maximum +velocity. + + courant = .5d0 + +[bool]normaladvect is a flag used to determine which algorithm to use to +calculate the new geoletry of the normals to the surfaces at the nodes on the +surfaces +if normaladvect is T, the normals are advected using the velocity gradient +if normaladvect is F, the normals are re-computed from the geometry of the +surface + + normaladvect = T + + +GRID ITERATIONS +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]griditer is a flag that allows for nonlinear iterations; when positive, a +fixed number (griditer) of iterations is permitted; when negative, the number of +nonlinear iterations is determined by a convergence criterion. + + griditer = -100 + +[dp]octree_refine_ratio is the threshold value used to determine whether the +octree has converged or not. the larger the value, the less stringent the test. + + octree_refine_ratio = .025d0 + + +NONLINEAR ITERATIONS +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]nonlinear_iterations is the maximum number of nonlinear iterations (i.e. +the iterations on a given constant grid) +if nonlinear_iterations is positive, it simply is the number of nonlinear +iterations performed for each grid. When negative it indicates an upper bound +of nonlinear iterations, but the actual number of nonlinear iterations is +determined by a convergence criterion (see the 'tol' parameter) + + nonlinear_iterations = -100 + +[int]nb_iter_nl_min is the number of nonlinear iterations required for velocity +convergence when using the convergence criterion tolerance value below. Note +that this is different than simply requiring a given number of nonlinear +iterations as the velocity field must also be converged within the tolerance +below. + + nb_iter_nl_min = 0 + +[dp]tol is the relative tolerance used to estimate convergence on the computed +velocity field + + tol = 0.01d0 + +[bool]adaptive_tol is a flag that toggles on/off the evolution of the tol +parameter with the grid level: when velocity convergence is reached on a grid, +the latest meaure of the velocity difference between the two last obtained +solutions is put in tol, so that on the following generated grid, the solution +reaches at least the same level of convergence. It allows to start with a not +too stringent value of tol at uniform octree level that evolves with every +grid, assuming that increasing the level of refinement of the octree allows to +better capture the solution, hence allowing a tighter convergence. + + adaptive_tol = F + + +OCTREES +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]leveluniform_oct is the level of uniform discretization of space; note that +a level is a power of two used to divide the unit cube + + leveluniform_oct = 4 + +[int]levelmax_oct is maximum level of octree discretization + + levelmax_oct = 5 + +[double]vex is the vertical exaggeration scaling factor that allows for +variable aspect ratio elements. The value should range between 0-1, which +will be multiplied by the vertical scale of the model domain. For example, a +value of 0.25 would correspond to a model of dimensions 1x1x0.25, where the +elements are 1/4 as tall as they are wide. NOTE: vex MUST be larger than the +highest surface/element in the model at all times. + + vex = 1.d0 + +[bool]ismooth is a flag to impose an additional level of smoothing after +refinement for the surfaces and strain rate. It ensures that no leaf is flanked +by other leaves differing by more than 1 level of refinement + + ismooth = F + +[int]noctreemax is the maximum size of any octree used in all computations + + noctreemax = 10000000 + +[dp]refine_ratio is used to determine octree refinement based on a given +criterion. All leaves where the criterion is larger than refine_ratio times the +maximum of this criterion are refined + + refine_ratio = -5500 + +[int]refine_criterion determines which refinement algorithm is to be used. +Several criteria exist for the refinement of the osolve octree. 1 is the second +invariant of the deviatoric strain-rate tensor; 2 is the sum of the squares of +the diagonal terms of the deviatoric strain-rate tensor; 3 is the second +invariant of the deviatoric strain rate tensor timses the leaf size. any other +value sets the criterion to zero and leads to no refinement. + + refine_criterion = 0 + +[int]initial_refine_level is the initial level at which the refinement of the +octree will be performed. it has to be smaller than levelmax_oct +this should be used (different from levelmax_oct) in case the flow is very +localized (nonlinear/plastic analysis) + + initial_refine_level = 4 + +[bool]renumber_nodes is a flag that can toggle on/off the renumbering of the +nodes by mean of Sloan's algorithm (T/F) + + renumber_nodes = T + + +BOUNDARY CONDITIONS AND VELOCITY CONSTRAINTS (much more coming in v0.2) +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[bool]2d_invariants forces the calculation of the strain rate tensor invariants +to two dimensions. This assumes no strain in the y direction, and sets eyy, exy +and eyz to zero. + + invariants_2d = F + + +PRESSURE +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]smoothing_type is a parameter allows to choose which type of smoothing is +to be applied to the pressure field: 0 is none, 1 is center->nodes->center, 2 is +the same, but weighted by neighbouring elemental volumes, 3 is regular grid+SPH, +and 4 is SPH. + + smoothing_type = 1 + + +CLOUD +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]npmin and [int]npmax are used to update the 3D volumetric cloud. npmin +corresponds to the minimum number of particles in any leaf; npmax is the maximum +allowable number in any leaf + + npmin = 1 + npmax = 4 + + +FEM + DIVFEM + MUMPS +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]matrule determines the rule used to assign material properties to +elements that are cut by a surface. The default value (0) will use divFEM to +determine the element properties. A value of 1 will assign the properties of +the material that occupies the majority of the element, whereas a value of 2 +will assign the properties of the volumetric minority of the element. NOTE: +divFEM is always used for any element containing void material. + + matrule = 1 + +[int]levelcut is the number of levels used to estimate the volume integrals in +the divFEM algorithm for cut cells; by testing, it has been estimated that a +value of 2 is very accurate + + levelcut = 2 + +[int]levelapprox is the number of levels used to estimate the remaining +integrals using an improved version of Marthijn's clever algorithm, usually 3 is +plenty + + levelapprox = 3 + +[dp]penalty is a global penalty parameter used to impose the bad faces or +incompatible faces linear constraints + + penalty = 1.d8 + +[bool]excl_vol is a parameter that toggles off the assumption that lsf's are +built on top of one another + + excl_vol = F + + +TEMPERATURE +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[dp] ztemp is the height interval between which a linear temperature gradient is +set: the temperature is 1 at the bottom, and 0 at ztemp. +[dp] tempscale + + calculate_temp = F + + ztemp = 0.2d0 + + tempscale = 100.d0 + + +ISOSTASY +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- +[bool]isostasy enables calculation of local isostasy with the equilibrium state +based on the initial geometry. Local isostasy will be used if this is true and +flexure is false. + + isostasy = F + +[bool]flexure enables calculation of flexural isostasy using the parameters +below. Note that either local or flexural isostasy will be used, based on this +value. Both must be true for flexure to be used. + + flexure = F + +[dp]elastic_plate_thickness is the thickness of the elastic plate used in the +flexural isostasy calculation, in meters (?) + + elastic_plate_thickness = 20.d3 + +[dp]density_difference is the difference in density between the crust and mantle +in ??? I believe this is abs(rhoc-rhom)*g. (kg/m^3 * m/s^2) + + density_difference = 3.d3 + +[bool]isobc enables modification of the velocity boundary conditions for certain +boundary condition types to account for changes in the Moho depth + + isobc = F + +MATERIALS +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]nmat is number of materials + + nmat = 1 + +[int]material0 determines what is the material above the 1st surface (the free +surface) if material0 is 0 then it is the void and the properties "0" are used +for that part of the model; otherwise the material is one of the material, +comprised between 1 and nmat + + material0 = 0 + +[bool]bulkvisc determines whether the user will input a material penalty value +that scales with the shear viscosity or whether an independent bulk viscosity +will be listed. If true, the values listed for [dp]penaltyi below should be +bulk viscosity values for each material. Otherwise, the listed penalty value +will be multiplied by the shear viscosity to enforce the incompressibility +condition. + + bulkvisc = F + +[bool]init_e2d determines whether an initial e2d value will be used to put +plastic materials on yield in the first time step, or whether the standard +initial material viscosity will be used. If true, an initial e2d value will be +used to ensure plastic materials behave as desired. + + init_e2d = T + +[dp]e2d0 is the initial strain rate value to use if init_e2d is true. This +value will be used in only the first nonlinear iteration of the first time +step and should be in scaled DOUAR units. + + e2d0 = 1.d-2 + +[dp]densityi, [dp]viscosityi and [dp]penaltyi are the density, viscosity and +incompressibility used for material i; there should be nmat sets of material +properties; there should also be a nil material if material0 has been set to 0 +[dp]expon is the nonlinear viscosity exponent +[dp]diffusivity is the heat diffusivity +[dp]expansion is the coefficient of thermal expansion +[dp]heat is the heat production +[char]plasticity_type is the type of plasticity +- No -> no plasticity, purely (nonlinear) viscous material +- vM -> von Mises yield criterion + -> [dp]plasticity_1st_param is the yield strength +- DPI-VII-> Drucker-Prager type of yield criterion + DPI -> the yield locus passes through the outer apices of the Mohr-Coulomb + hexagonal cone + DPII -> the yield locus passes through the inner apices of the Mohr-Coulomb + hexagonal cone + DPIII -> The average case, intermediate to DPI and DPII (USE ONLY FOR TESTING; + this case can not easily be derived as a simplification of the + Mohr-Coulomb criterion) + DPIV -> Plane strain Drucker-Prager formulation for associative materials + (USE ONLY FOR TESTING) + -> [dp]plasticity_1st_param is the angle phi + -> [dp]plasticity_2nd_param is the cohesion c + -> [dp]plasticity_3rd_param is the accumulated strain minimum, beyond + which phi is modified + -> [dp]plasticity_4th_param is the accumulated strain maximum, beyond + which phi is set to the value in plasticity_5th_param + -> [dp]plasticity_5th_param is the final phi value, if modified +- DPV -> Drucker-Prager type of yield criterion that does not calculate the + values of alpha and k from phi and c. They are specified here. Note + that strain softening cannot be used for this DP type. + -> [dp]plasticity_1st_param is alpha + -> [dp]plasticity_2nd_param is k +- DPVI -> Drucker-Prager type of yield criterion that uses the equation below + sigma_y = pressure * tan(phi) + c + (USE ONLY FOR TESTING) + -> [dp]plasticity_1st_param is the angle phi + -> [dp]plasticity_2nd_param is the cohesion c + -> [dp]plasticity_3rd_param is the accumulated strain minimum, beyond + which phi is modified + -> [dp]plasticity_4th_param is the accumulated strain maximum, beyond + which phi is set to the value in plasticity_5th_param + -> [dp]plasticity_5th_param is the final phi value, if modified +- DPVII -> Drucker-Prager type of yield criterion that uses the same 2D + formulation that is used in Sopale + sigma_y = pressure*sin(phi) + c*cos(phi) + -> [dp]plasticity_1st_param is the angle phi + -> [dp]plasticity_2nd_param is the cohesion c + -> [dp]plasticity_3rd_param is the accumulated strain minimum, beyond + which phi is modified + -> [dp]plasticity_4th_param is the accumulated strain maximum, beyond + which phi is set to the value in plasticity_5th_param + -> [dp]plasticity_5th_param is the final phi value, if modified +- MC -> Mohr-Coulomb type of yield criterion + -> [dp]plasticity_1st_param is the angle phi + -> [dp]plasticity_2nd_param is the cohesion c + -> [dp]plasticity_3rd_param is the accumulated strain minimum, beyond + which phi is modified + -> [dp]plasticity_4th_param is the accumulated strain maximum, beyond + which phi is set to the value in plasticity_5th_param + -> [dp]plasticity_5th_param is the final phi value, if modified +[dp]ztrans is the elevation below which the material transition from the given +material to material transnum +[int]transnum is the material number used when z is below ztrans + + density0 = 0.d0 + viscosity0 = 1.d-5 + penalty0 = 1.d8 + expon0 = 1.d0 + diffusivity0 = 1.d0 + heat0 = 0.d0 + activationenergy0 = 0.d0 + plasticity_type0 = No + plasticity_1st_param0 = 0.d0 + plasticity_2nd_param0 = 0.d0 + plasticity_3rd_param0 = 0.5d0 + plasticity_4th_param0 = 1.5d0 + plasticity_5th_param0 = 0.d0 + ztrans0 = -1.d0 + transnum0 = 0 + + density1 = -1.d0 + viscosity1 = 5.d4 + penalty1 = 1.d8 + expon1 = 1.d0 + diffusivity1 = 2.d-3 + heat1 = 0.d0 + activationenergy1 = 0.d0 + plasticity_type1 = MC + plasticity_1st_param1 = 15.d0 + plasticity_2nd_param1 = 5.d-4 + plasticity_3rd_param1 = 0.5d0 + plasticity_4th_param1 = 1.5d0 + plasticity_5th_param1 = 15.d0 + ztrans1 = -1.d0 + transnum1 = 1 + +[dp]viscositymin and viscositymax are bounds on the viscosity (if negative bound +is not imposed) These bounds are introduced to prevent the viscosity to reach +unrealistic values, especially when using non-linear (power-law or brittle) +rheologies + + viscositymin = 5.d-8 + viscositymax = 5.d9 + + +SURFACES +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]ns is number of surfaces to track + + ns = 1 + +for each surface, one needs to define a levelt, itype, material, surface_type, +activation_time and surface_params. +- [int]levelt is the inital level for the particles on the surface; to be + accurate and avoid holes in the surface during definition of the lsf, one + should use levelt eq to levelmax_oct+1 for all surfaces as a minimum value; +- [int]itype should be 1 for foldable surfaces or 0 for nonfoldable surfaces; +- [int]material is the material type refering to the table of material available + (max nmat); +- [dp]activation_time is the time the surface becomes active (before that time, + it is glued to the 0th surface). this parameter is useful when defining + stratigraphic horizons; default is -1, ie the surface is not glued to the free + surface +- [int]surface_type is comprised between 1 and 13. + 1 corresponds to a flat surface, + -> surface_param_01 is the z level + 2 to a rectangular emboss, + -> surface_param_01 is the z level + -> surface_param_02 and 03 are x1,x2 + -> surface_param_04 and 05 are y1,y2 + -> surface_param_06 is the thickness + 3 to a convex spherical emboss, + -> surface_param_01 is the z level + -> surface_param_02 and 03 are x0,y0 + -> surface_param_04 is the radius + 4 to concave spherical emboss, + -> surface_param_01 is the z level + -> surface_param_02 and 03 are x0,y0 + -> surface_param_04 is the radius + 5 to a double rectangular emboss, + -> surface_param_01 is the z level + -> surface_param_02 and 03 are x1,x2 + -> surface_param_04 and 05 are x3,x4 + -> surface_param_06 and 07 are y1,y2 + -> surface_param_08 and 09 are y3,y4 + -> surface_param_10 is the thickness + 6 to a sinus, + -> surface_param_01 is the z level + -> surface_param_02 is the wavelength + -> surface_param_03 is the amplitude + 7 to a noisy surface, + -> surface_param_01 is the z level + -> surface_param_02 is the noise amplitude + 8 to a double sinus, + -> surface_param_01 is the z level + -> surface_param_02 is the x-wavelength + -> surface_param_03 is the x-amplitude + -> surface_param_04 is the y-wavelength + -> surface_param_05 is the y-amplitude + 9 to a cosinus, + -> surface_param_01 is the z level + -> surface_param_02 is the wavelength + -> surface_param_03 is the amplitude + 10 to a 2D embankment, + -> surface_param_01 is z0 + -> surface_param_02 is y0 + -> surface_param_03 is psi + -> surface_param_04 is the thickness + 11 to a 2D hill, + -> surface_param_01 is z level + -> surface_param_02 is the height + -> surface_param_03 is y0 + -> surface_param_04 is width + 12 to a rectangular emboss with specified slope, + -> surface_param_01 is the z level + -> surface_param_02 is x1 + -> surface_param_03 is x2 + -> surface_param_04 is y1 + -> surface_param_05 is y2 + -> surface_param_06 is the thickness + -> surface_param_07 is the slope + 13 to a 2D embankment of finite length and with 2 kinks + -> surface_param_01 is the z level + -> surface_param_02 is the thickness (zthick) + -> surface_param_03 is slope angle (psi) + -> surface_param_04 is x value for the center of the slope (x0) + -> surface_param_05 is the x width of the narrow end of the plateau (y0w) + -> surface_param_06 is the y position for the start of the slope (y0) + -> surface_param_07 is the y position of the first kink (y1) + -> surface_param_08 is the angle of the first kink (theta1) + -> surface_param_09 is the y position of the second kink (y2) + -> surface_param_10 is the angle of the second kink (theta2) + +- [int]leveloct is the level at which the octree will be refined in the vicinity + of the surface. + + levelt1 = 6 + itype1 = 1 + surface_type_1 = 1 + rand1 = F + surface_param_01_1 = 0.2d0 + material1 = 1 + activation_time_1 = -1. + leveloct1 = 5 + stretch1 = 1.5d0 + anglemax1 = 180.d0 + criterion1 = 1 + anglemaxoctree1 = 180.d0 + spread_surface_points1 = 1 + +[int]niter_move is the number of iterations used to update particle positions +using an implicit, mid-point algorithm (default is 10) + + niter_move = 10 + +[dp]stretch is the maximum allowed increase in linear length between two +initially adjacent particles on any surface; when this stretch is achieved, a +new particle is inserted on the surface, half-way along the stretched edge + + stretch = 1.5d0 + +[dp]anglemax is the maximum allowed angle between two adjacent normals +when the angle is reached a new point is inserted bewteen the two points to +reduce the angle between the two normals + + anglemax = 180.d0 + +[int]criterion is criterion used to define the octree in the vicinity of the +sufaces; criterion 1 corresponds to imposing that all leaves that are cut by any +of the surfaces must be at level levelmax_oct; criterion 2 corresponds to +imposing that discretization is proportional to the curvature of the surface; +curvature is calculated from the local divergence of the normals. criterion 3 +corresponds to imposing that all leaves that contain at least one particle of +any surface is at levelmax_oct; + + criterion = 2 + +[dp]anglemaxoctree is only defined for criterion 2; t is the maximum allowable +angle between two adjacent normals; if the angle is greater than anglemaxoctree, +the local octree leaves are forced to be at level levelmax_oct; otherwise they +are proportionally larger (smaller levels) (default is 10) + + anglemaxoctree = 10.d0 + + +REFINEMENT IN BOXES +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]nboxes is the number of boxes in which the user imposes a set level of +discretization + + nboxes = 1 + +for each box we need to specify the two end corners and the level +the syntax is [dp]boxnx0, [dp]boxnx1, [dp]boxny0, [dp]boxny1, [dp]boxnz0, +[dp]boxnz1 and boxnlevel where n is the box number + + box1x0 = 0.d0 + box1x1 = 1.d0 + box1y0 = 0.d0 + box1y1 = 1.d0 + box1z0 = 0.d0 + box1z1 = 0.25d0 + box1level = 5 + + +REFINEMENT ON CUBE FACES +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- +[bool]ref_on_faces toggles on/off the user imposed refinement on faces. +For each of the six faces of the cube, on defines the level at which the +desired area is to be refined. this area is given by bottom, top, left and +right coordinates that are comprised between 0 and 1 +faces 1,2,3,4,5,6 respectively correspond to x=0,x=1,y=0,y=1,z=0,z=1 + + ref_on_faces = F + + level_face1 = 5 + b1 = .02 + t1 = .51 + l1 = .11 + r1 = .81 + level_face2 = 5 + b2 = .45 + t2 = .55 + l2 = .46 + r2 = .56 + level_face3 = 5 + b3 = .0 + t3 = .2 + l3 = .3 + r3 = .7 + level_face4 = 5 + b4 = .4 + t4 = .5 + l4 = .14 + r4 = .4 + level_face5 = 6 + b5 = 0. + t5 = 1. + l5 = 0. + r5 = 0.75 + level_face6 = 5 + b6 = .26 + t6 = .56 + l6 = .16 + r6 = .86 + + +EROSION +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[bool]erosion is a flag that toggles on/off the erosion. + + erosion = T + +if erosion is on, one also needs to set the erosion level/height, [dp]zerosion +Note that this is a first attempt at erosion; in future versions, DOUAR +should be easily linked to a surface processes model like CASCADE + + zerosion = 0.2d0 + +If erosion is on, one also needs to define a length scale and a velocity scale +to properly translate the uplift rate produced by DOUAR into something that is +adequate for CASCADE; +[dp]length_scale is the scale of the unit model in km +[dp]velocity_scale is the velocity scale in km/Myr. +if length_scale is negative, erosion is assumed to be perfect (no call to +CASCADE is needed) + + length_scale = 640.d0 + velocity_scale = 50.d0 + +One then needs to define the erosion constants in CASCADE +[dp]fluvial_erosion is the fluvial erosion constant in 1/m^2 +[dp]diffusion_erosion is the diffusion erosion constant in m^2/yr +(4d-2 32d-2) + + fluvial_erosion = 4.d-2 + diffusion_erosion = 32.d-2 + +One finally needs to specify the boundary conditions for CASCADE +if [int]baselevelx0 is set to 1 then the boundary at x=0 is set at baselevel +(water and sediment exit) +if [int]baselevelx1 is set to 1 then the boundary at x=1 is set at baselevel +(water and sediment exit) +if [int]baselevely0 is set to 1 then the boundary at y=0 is set at baselevel +(water and sediment exit) +if [int]baselevely1 is set to 1 then the boundary at y=1 is set at baselevel +(water and sediment exit) + + baselevelx0 = 1 + baselevelx1 = 1 + baselevely0 = 0 + baselevely1 = 1 + + +MATRIX VISUALISATION +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[bool]visualise_matrix allows the user to turn on the visual representation of +the matrices used in the code. Be careful, the generated postscript files are +huge! + + visualise_matrix = F + + +CROSS SECTIONS +_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_- + +[int]nsections is the number of cross-sections to be output. [dp]scale is the +scale used to produce the postscript file. +[int]xyz takes values 1, 2 or 3, and corresponds to planes defined by +x=constant, y=constant, and z=constant resp. +The [bool]flags are self explanatory. [char]colormap is the chosen colormap (jet +or hot). [int]ncolours is the number of colours used to produce the plot. + +nsections = 0 + +xyz_1 = 1 +slice_1 = 0.5001 +flag_press_1 = T +flag_spress_1 = T +flag_e2d_1 = T +flag_e3d_1 = F +flag_strain_1 = F +flag_lode_1 = F +flag_crit_1 = F +flag_grid_1 = T +flag_mu_1 = F +flag_u_1 = F +flag_v_1 = F +flag_w_1 = F +flag_q_1 = F +flag_uvw_1 = F +flag_lsf_1 = F +flag_vfield_1 = F +flag_colour_1 = T +flag_plastic_1 = F +flag_temp_1 = F +flag_velvect_1 = T +scale_1 = 500. +colormap_1 = jet +ncolours_1 = 256 + +xyz_2 = 2 +slice_2 = 0.901 +flag_press_2 = T +flag_e2d_2 = T +flag_e3d_2 = F +flag_strain_1 = F +flag_lode_2 = F +flag_crit_2 = F +flag_grid_2 = F +flag_mu_2 = F +flag_u_2 = F +flag_v_2 = F +flag_w_2 = F +flag_q_2 = F +flag_uvw_2 = F +flag_lsf_2 = F +flag_vfield_2 = F +flag_colour_2 = T +flag_plastic_2 = F +flag_temp_1 = F +flag_velvect_2 = F +scale_2 = 800. +colormap_2 = jet +ncolours_2 = 256 + +xyz_3 = 3 +slice_3 = 0.0010 +flag_press_3 = F +flag_e2d_3 = F +flag_e3d_3 = F +flag_strain_1 = F +flag_lode_3 = F +flag_crit_3 = F +flag_grid_3 = F +flag_mu_3 = F +flag_u_3 = F +flag_v_3 = F +flag_w_3 = F +flag_q_3 = F +flag_uvw_3 = F +flag_lsf_3 = F +flag_vfield_3 = F +flag_colour_3 = T +flag_plastic_3 = F +flag_temp_1 = F +flag_velvect_3 = F +scale_3 = 400. +colormap_3 = jet +ncolours_3 = 256