source: trunk/anuga_core/source/anuga/config.py @ 7866

Last change on this file since 7866 was 7865, checked in by hudson, 14 years ago

Refactoring to clean up pylint errors.

File size: 9.3 KB
Line 
1"""Module where global ANUGA model parameters and default values are set
2"""
3
4import os
5import sys
6
7
8################################################################################
9# Numerical constants
10################################################################################
11
12epsilon = 1.0e-12                    # Smallest number - used for safe division
13max_float = 1.0e36                   # Largest number - used to initialise
14                                     # (max, min) ranges
15default_smoothing_parameter = 0.001  # Default alpha for penalised
16                                     # least squares fitting
17single_precision = 1.0e-6            # Smallest single precision number
18velocity_protection = 1.0e-6         # Used to compute velocity from momentum
19                                     # See section 7.4 on Flux limiting
20                                     # in the user manual
21                           
22
23################################################################################
24# Standard filenames, directories and system parameters used by ANUGA
25################################################################################
26
27pmesh_filename = '.\\pmesh'
28version_filename = 'stored_version_info.py'
29default_datadir = '.'
30time_format = '%d/%m/%y %H:%M:%S'    # Used with timefile2netcdf
31umask = 002  # Controls file and directory permission created by anuga (UNIX)
32default_boundary_tag = 'exterior' 
33
34# Major revision number for use with create_distribution
35# and update_anuga_user_guide
36major_revision = '1.2'
37
38################################################################################
39# Physical constants
40################################################################################
41
42manning = 0.03  # Manning's friction coefficient
43#g = 9.80665    # Gravity - FIXME reinstate this and fix unit tests.
44g = 9.8
45#g(phi) = 9780313 * (1 + 0.0053024 sin(phi)**2 - 0.000 0059 sin(2*phi)**2)
46# micro m/s**2, where phi is the latitude
47#The 'official' average is 9.80665
48
49eta_w = 3.0e-3 # Wind stress coefficient
50rho_a = 1.2e-3 # Atmospheric density
51rho_w = 1023   # Fluid density [kg/m^3] (rho_w = 1023 for salt water)
52
53################################################################################
54# Limiters - used with linear reconstruction of vertex
55# values from centroid values
56################################################################################
57
58# Betas [0;1] control the allowed steepness of gradient for second order
59# extrapolations. Values of 1 allow the steepes gradients while
60# lower values are more conservative. Values of 0 correspond to
61# 1'st order extrapolations.
62#
63
64# There are separate betas for the w, uh, and vh limiters
65# I think these are better SR but they conflict with the unit tests!
66beta_w      = 1.0
67beta_w_dry  = 0.2
68beta_uh     = 1.0
69beta_uh_dry = 0.2
70beta_vh     = 1.0
71beta_vh_dry = 0.2
72
73# Alpha_balance controls how limiters are balanced between deep and shallow.
74# A large value will favour the deep water limiters, allowing the a closer hug
75# to the coastline.  This will minimise 'creep' but at the same time cause
76# smaller time steps
77# Range:
78alpha_balance = 2.0 
79
80# Flag use of new limiters.
81# tight_slope_limiters = 0 means use old limiters (e.g. for some tests)
82# tight_slope_limiters = 1 means use new limiters that hug the bathymetry closer
83tight_slope_limiters = True
84
85# Use centroid velocities to reconstruct momentum at vertices in
86# very shallow water
87# This option has a first order flavour to it, but we still have second order
88# reconstruction of stage and this option only applies in
89# balance_deep_and_shallow when
90# alpha < 1 so in deeper water the full second order scheme is used.
91#
92# This option is good with tight_slope_limiters, especially for large domains.
93use_centroid_velocities = True
94       
95# FIXME (Ole) Maybe get rid of order altogether and use beta_w
96default_order = 2
97
98################################################################################
99# Timestepping
100################################################################################
101
102CFL = 1.0  # CFL condition assigned to domain.CFL - controls timestep size
103     
104# Choose type of timestepping,
105#timestepping_method = 'rk2'   # 2nd Order TVD scheme
106timestepping_method = 'euler' # 1st order euler
107
108# rk2 is a little more stable than euler, so rk2 timestepping
109# can deal with a larger beta when slope limiting the reconstructed
110# solution. The large beta is needed if solving problems sensitive
111# to numerical diffusion, like a small forced wave in an ocean
112beta_euler = 1.0
113beta_rk2   = 1.6
114
115# Option to search for signatures where isolated triangles are
116# responsible for a small global timestep.
117# Treating these by limiting their momenta may help speed up the
118# overall computation.
119# This facility is experimental.
120# protect_against_isolated_degenerate_timesteps = False
121protect_against_isolated_degenerate_timesteps = False
122
123min_timestep = 1.0e-6 # Minimal timestep accepted in ANUGA
124max_timestep = 1.0e+3
125max_smallsteps = 50   # Max number of degenerate steps allowed b4
126                      # trying first order
127
128# Perhaps minimal timestep could be based on the geometry as follows:
129# Define maximal possible speed in open water v_max, e.g. 500m/s (soundspeed?)
130# Then work out minimal internal distance in mesh r_min and set
131# min_timestep = r_min/v_max
132#
133# Max speeds are calculated in the flux function as
134#
135# lambda = v +/- sqrt(gh)
136#
137# so with 500 m/s, h ~ 500^2/g = 2500 m well out of the domain of the
138# shallow water wave equation
139#
140# The actual soundspeed can be as high as 1530m/s
141# (see http://staff.washington.edu/aganse/public.projects/
142#            clustering/clustering.html),
143# but that would only happen with h>225000m in this equation. Why ?
144# The maximal speed we specify is really related to the max speed
145# of surface pertubation
146#
147# v_max = 100 #For use in domain_ext.c
148# sound_speed = 500
149
150################################################################################
151# Ranges specific to the shallow water wave equation
152# These control maximal and minimal values of quantities
153################################################################################
154
155# Water depth below which it is considered to be 0 in the model
156minimum_allowed_height = 1.0e-3 
157
158# Water depth below which it is *stored* as 0
159minimum_storable_height = 1.0e-5
160
161# FIXME (Ole): Redefine this parameter to control maximal speeds in general
162# and associate it with protect_against_isolated_degenerate_timesteps = True
163maximum_allowed_speed = 0.0 # Maximal particle speed of water
164#maximum_allowed_speed = 1.0 # Maximal particle speed of water
165                             # Too large (100) creates 'flopping' water
166                             # Too small (0) creates 'creep'
167
168maximum_froude_number = 100.0 # To be used in limiters.
169
170################################################################################
171# Performance parameters used to invoke various optimisations
172################################################################################
173
174use_extensions = True # Use C-extensions
175use_psyco = True      # Use psyco optimisations
176
177optimise_dry_cells = True # Exclude dry and still cells from flux computation
178optimised_gradient_limiter = True # Use hardwired gradient limiter
179use_edge_limiter = False  # The edge limiter is better, but most runs have been
180                          # using vertex limiting. Validations passed with this
181                          # one True 9th May 2008, but many unit tests need
182                          # backward compatibility flag set FIXME(Ole).
183
184points_file_block_line_size = 500 # Number of lines read in from a points file
185                                  # when blocking
186
187################################################################################
188# NetCDF-specific type constants.  Used when defining NetCDF file variables.
189################################################################################
190
191netcdf_char = 'c'
192netcdf_byte = 'b'
193netcdf_int = 'i'
194netcdf_float = 'd'
195netcdf_float64 = 'd'
196netcdf_float32 = 'f'
197
198################################################################################
199# Dynamically-defined constants.
200################################################################################
201
202# Determine if we can read/write large NetCDF files
203netcdf_mode_w = 'w'
204netcdf_mode_a = 'a'
205netcdf_mode_r = 'r'
206
207# Code to set the write mode depending on
208# whether Scientific.IO supports large NetCDF files
209s = """
210import os, tempfile
211from Scientific.IO.NetCDF import NetCDFFile
212
213filename = tempfile.mktemp('.nc')
214
215fid = NetCDFFile(filename, 'wl')
216fid.close()
217os.remove(filename)
218"""
219
220# Need to run in a separate process due an
221# error with older versions of Scientific.IO
222if sys.platform == 'win32':
223    null = 'NUL'
224else:
225    null = '/dev/null'
226cmd = 'python -c "%s" 2> %s' % (s, null)
227err = os.system(cmd)
228
229if err != 0:
230    # The Python script s failed e.g. with a segfault
231    # which means that large file support is
232    # definitely not supported
233    pass
234else:   
235    # Try the import within this process
236    try:
237        exec(s)
238    except IOError:
239        # NetCDFFile does not segfault but it does not
240        # support large file support   
241        pass
242    else:
243        # Set the default mode to large file support
244        #netcdf_mode_w = 'w4' # Future use of HDF5       
245        netcdf_mode_w = 'wl' # Large NetCDF (new default 30/6/2009)
246        #netcdf_mode_w = 'w'   # Old style NetCDF used by OSG viewer
247
248
Note: See TracBrowser for help on using the repository browser.