1. 01 Nov, 2017 1 commit
    • Christopher Subich's avatar
      Make timing code compile-time optional (TIMINGS=true) · d68f6d2b
      Christopher Subich authored
      This changeset adds a compile-time option on whether to include the
      timing code.  If SPINS is compiled without TIMINGS=true, the timing code
      is #defined to a no-op and should not otherwise affect the performance
      of the software; timing_stack_report will instead print a short error
      message.
      
      Because this affects #defines in header files, changing from
      TIMINGS=true to false (or vice versa) should only be done with a 'make
      clean'.  Not doing so could leave the timing code in an inconsistent
      state and result in either link-time errors (if older code calls timing
      code that is #defined away) or a corrupted/growing timing stack (if for
      example a timing_push() is compiled in but the corresponding
      timing_pop() is removed.)
      d68f6d2b
  2. 26 Oct, 2017 1 commit
    • Christopher Subich's avatar
      Simplified tridiangonal solve algorithm; timing extensions · 3581d1c1
      Christopher Subich authored
      First and most trivially, this change adds a few more timing push/pop
      invocations relating to multigrid and especially to some of the MPI
      synchronization steps that might cause processors to wait on each other.
      
      More substantively, this change also adds a new tridiangonal solver
      based on the bog-standard forwards/backwards elimination Thompson
      algorithm (see Wikipedia).  This should be acceptable because the 1D
      problems being solved themselves come from a 2D problem, so we don't
      expect ill-conditioning; calling out to the GMRES banded solver was
      surprisingly a computational bottleneck perhaps because of pivoting.
      
      This change seems to decrease the line-solve time by about 80%, which in
      turn decreases the overall runtime (tank_rho test case) by 40%.
      3581d1c1
  3. 29 Jun, 2017 2 commits
  4. 26 Jun, 2017 1 commit
    • Christopher Subich's avatar
      Stack-structured timing code · 0a4781c9
      Christopher Subich authored
      This commit adds stack-structured timing code, such that bits of SPINS
      can be instrumented with
      
      >  timing_push("short name")
      >  [... work goes here]
      >  timing_pop()
      
      This is recursive, such that already-timed portions of the code can call
      the timing code to further instrument subroutines or discrete code
      segments, and the resulting structure follows the -dynamic- call stack.
      This does have the odd result that portions of the code that are re-used
      (namely spectral derivatives) can show up as several 'leaf' nodes in the
      graph.
      
      To print out the report to standard error, call
      
      >  timing_stack_report()
      
      The exact format of the report should be considered tentative.
      0a4781c9
  5. 15 Jun, 2017 1 commit
  6. 25 Apr, 2017 1 commit
    • David Deepwell's avatar
      Add case for computing derivatives and vorticity · 92b79931
      David Deepwell authored
      A new case file, derivatives.cpp, gives a means to compute
      derivatives of any field or to compute vorticity components.
      The vorticity calculation in Science is completely re-written
      so-as to also work for mapped grids.
      
      The derivatives file is built to compute secondary variables
      after a run has already been completed. This uses spins' built-in
      derivative toolkit allowing multiple variables to be calculated
      in parallel. A derivative field (one that is the derivative of
      another) is denoted by an underscore followed by the direction
      the derivative was taken in (ie. u_x is the x derivative of u).
      92b79931
  7. 22 Mar, 2017 1 commit
    • Christopher Subich's avatar
      Make default get_dt_max error on invocation · 075caad6
      Christopher Subich authored
      The method get_dt_max inside BaseCase was added in a
      previous commit as part of the effort to include the
      'default' timestep-checking code inside BaseCase.
      
      The intent is for user code to override this method with a
      problem-specific maximum timestep, but if the user code did
      not do this the previous implementation would return '0'.
      
      This would cause SPINS to error, as a timestep of <= 0 is
      obviously invalid.  However, the error message would be
      unintuivie, referring to the <=0 timestep rather than the
      root cause of the user not implementing get_dt_max().
      
      This change replaces the previous default method with one
      that immediately aborts, with an assertion failure if
      assert() is enabled at compile time.
      075caad6
  8. 04 Feb, 2017 3 commits
    • David Deepwell's avatar
      Move time step checker to BaseCase · 428b5399
      David Deepwell authored
      The empty shell of check_timestep in BaseCase has been filled
      with what was essentially the same thing in each case file.
      Again, this is just more case file house cleaning.
      
      A few ancillary functions have also been written to help.
      The most important might be get_dt_max() which returns dt_max.
      dt_max can either be specified in spins.conf or defined in the
      case file based on the buoyancy frequency. Currently, if
      dt_max > 0 and in spins.conf, then that value is used. Otherwise,
      0.5/sqrt(N2_max) is used.
      
      The re-definition of check_timestep in each case file can now be
      removed so long as the associated ancillary functions and
      parameters are included.
      428b5399
    • David Deepwell's avatar
      Move restart sequence checking into Options.cpp · 25b5d62f
      David Deepwell authored
      The generic code block to check the restart sequence is not needed
      in the case file and has been moved into Options.cpp. This
      streamlines the case file and makes a single copy for all other
      case files to use.
      25b5d62f
    • David Deepwell's avatar
      Properly initialize get_restart_sequence() in case file · a1a18e74
      David Deepwell authored
      The re-definition of get_restart_sequence() was missing in
      gravity_current.cpp. This resulted in all restarts to read output
      zero (*.0). Restarting will now properly read the output number
      specified by restart_sequence.
      a1a18e74
  9. 02 Feb, 2017 1 commit
    • Christopher Subich's avatar
      Move BC-parsing helper to BaseCase · e3adc1f9
      Christopher Subich authored
      This commit moves the helper function that process
      expansion/boundary conditions ("FOURIER" -> PERIODIC, etc)
      to BaseCase.cpp, from Options.
      
      This change should be invisible to most case files, which
      habitually include BaseCase for helper functions as-is, but
      it allows Options to once again be independent of
      NSIntegrator.  This improves modularity for other uses of
      the code, such as a hypothetical SPINSoff (get it?  it's a
      pun) that would like to do some spectral voodoo with an
      option file but doesn't need the full Navier-Stokes
      machinery.
      e3adc1f9
  10. 27 Jan, 2017 2 commits
    • David Deepwell's avatar
      Add filter parameters to case options · b5470525
      David Deepwell authored
      gravity_current.cpp accepts filter parameters (f_cutoff, f_order,
      f_strength) as optional input arguments. Easy adjustment of the
      filter is useful when the flow become suddenly becomes too chaotic.
      b5470525
    • David Deepwell's avatar
      Move expansion parsing into Options.cpp · 3035b35a
      David Deepwell authored
      Parsing of expansion types in gravity_current.cpp moved to
      Options.cpp. This is a general procedure that other case files
      will and should make use of. It also simplifies the case file
      to the specific set-up of that case.
      3035b35a
  11. 26 Jan, 2017 1 commit
    • Christopher Subich's avatar
      Update maximum-precision float format to .17g · 46ca2c72
      Christopher Subich authored
      See http://stackoverflow.com/questions/16839658 -- although
      in most cases 16 decimal digits suffices for a double
      precision value, in a handful of instances 17 are necessary.
      At the command line, see the difference between:
      
      $ printf %.17g 1.0000000000000001
      
      and
      
      $ printf %.16g 1.0000000000000001
      
      Since SPINS does not provide bit-for-bit compatibility when
      restarting (because of the startup timestepping approach),
      this change is unlikely to make any significant difference
      to existing or future cases.
      46ca2c72
  12. 24 Jan, 2017 3 commits
    • David Deepwell's avatar
      Increase accuracy of dump time. · 0208ed34
      David Deepwell authored
      The accuracy of the time used in dump_time.txt (the dump time)
      just before writing dump fields was increased from 12 significant
      figures to 16. 12 was probably fine, but since a double is used
      we might as well save all digits of accuracy.
      0208ed34
    • David Deepwell's avatar
      Moved parameter adjustment for restarting from dump from case file into Options.cpp · 577a2a1d
      David Deepwell authored
      Moved the section of gravity_current.cpp regarding the adjustment
      of temporal values (restarting, initial_time, etc.) when restarting
      from dump into Options.cpp as it is a standard piece of code and
      should be in a centralized location.
      577a2a1d
    • Christopher Subich's avatar
      Modification of gravity_current to support subdir build · f73b2120
      Christopher Subich authored
      The gravity_current case was located in a subdirectory, but
      its #include directives had only one ".." in the path.  This
      updates those directive to properly refer to the grandparent
      subdirectory (src/), allowing the case to be built as:
      
      make -j# cases/gravity_current/gravity_current.x
      f73b2120
  13. 16 Jan, 2017 1 commit
    • David Deepwell's avatar
      Creation of a new case file and associated configuration file. · 4a8d478c
      David Deepwell authored
      This case file is restartable and will dump all variables at the end of the requested computation time (stated in spins.conf) when run on SCINET or SHARCNET. The spins.conf enables the quick set-up of a parameter sweep as parameters are set within it and not within the case file which would need to be re-compiled.
      This file will be used as a basis for the upcoming additions and changes to BaseCase.cpp and Science.cpp.
      4a8d478c
  14. 09 Sep, 2016 1 commit
    • Christopher Subich's avatar
      Add vector-valued options to the config file · 9e8fc9c6
      Christopher Subich authored
      This change modifeis Options.hpp to add the possibility of
      vector-valued options to the config file and command line.
      From the 'case' code, this option is used by (for example)
      
      std::vector<int> intvec;
      add_option("Name",&intvec,"Description);
      
      When run, options of the form --Name 1 2 3 4 (at the command
      line) or Name = 1 2 3 4 (in the config file) will be added
      to the vector 'intvec' as separate entries.
      
      Because add_option is a template function, this works
      identically for double (floating point)-valued options.  It
      has not been tested with string options.
      
      Since an empty vector is a perfectly cromulent and
      unambiguous possibility, vector-valued options do *not*
      permit the specification of defaults.
      9e8fc9c6
  15. 21 Jun, 2016 2 commits
  16. 11 Mar, 2016 1 commit
    • Christopher Subich's avatar
      Add caching for array-rebalance in multigrid · 3164ef96
      Christopher Subich authored
      The multigrid solver regularly "rebalances" arrays, by which
      it takes an input array and moves complete rows (i fixed, j
      variable) such that the same global arrays divides
      differently.  This is used on initial calling (to adjust the
      input problem specification), on coarsening, and finally to
      collapse "down" processors when near the coarsest level.
      
      However, this code was written naïvely: it built (and
      divided) the global array split at each invocation, using a
      number of MPI_Allgather calls on each processor.  Profiling
      (courtesy Kris Rowe) has found that this can be agonizingly
      slow; in at least one test case on one particular system it
      consumed over 50% of wallclock-time.
      
      This patch fixes the problem by allowing rebalance_array
      (now moved into the multigrid solver class) to cache these
      array divisions.  It requires the caller to specify one of a
      few categories, with the division being computed (with
      Allgathers) only on the first invocation per multigrid
      level and type.  The divisions currently are:
      
      FFine, UFine, FCoarse, and UCoarse
      
      where 'F' and 'U' refer to the error and solution terms
      respectively.  Confusingly, coefficients of the Laplacian
      itself (such as the uxx term) belong to the 'error' space in
      terms of call structure.
      
      This caching eliminates all but setup calls to Allgather.
      3164ef96
  17. 15 Dec, 2015 1 commit
  18. 28 Oct, 2015 2 commits
    • Christopher Subich's avatar
      Fix compilation of automatic_grid · 1ac02e05
      Christopher Subich authored
      The automatic_grid() helper function in BaseCase was not properly
      functional if called without supplied arguments for the 1D grids
      xx/yy/zz.  Now, if these arguments are null or not supplied, the
      function will do its grid generation using temporarily-allocated
      internal 1D arrays, deleted at the end of its execution.
      1ac02e05
    • Christopher Subich's avatar
      Correct sign error in automatic generation of grids · 5c13e0d9
      Christopher Subich authored
      Previously, the automatic_grid function was creating Chebyshev-type
      grids with points ordered as:
      
      {x,z} = L/2 + L/2*cos(pi*ii/(N-1))
      
      This is incorrect.  The expected ordering is actually the reverse, or
      
      {x,z} = L/2 - L/2*cos(...)
      5c13e0d9
  19. 18 Sep, 2015 2 commits
  20. 24 Aug, 2015 1 commit
  21. 10 Jun, 2015 1 commit
  22. 28 Nov, 2014 1 commit
    • Christopher Subich's avatar
      Fix for wavenumbers of null transform · a41e044e
      Christopher Subich authored
      Parformer.cpp -- previously, an all-NULL transform being
      asked for its wavenumbers would only return an array of size
      one, rather than one of appropriate length for the (real)
      "temporary" array.  This could conceivably have caused
      problems for the graceful use of 3D routines on
      2D-but-layered fields.
      a41e044e
  23. 29 Sep, 2014 2 commits
    • Christopher Subich's avatar
      Removed old "Transformer" code · 1a1c8814
      Christopher Subich authored
      The pair of Transformer.[ch]pp files contained old code built for SPINS
      before parallel transformation (over MPI) was fully implemented.  These
      classes were replaced wholesale with Parformer.[ch]pp (Parallel
      Transformer), but the older code remained in the repository.
      
      Removing these files should not affect any active code.  It will affect
      certain disued test cases (subject to changing #includes), but these
      cases themselves need more serious maintenance to update for the
      "gradient object" interface.
      1a1c8814
    • Michael's avatar
      Print descriptive names for dimtype and expansion · c8cdf0c1
      Michael authored
      At startup, print a string describing the dimension
      type and spectral expansion type, rather than the
      internal enum value.
      c8cdf0c1
  24. 18 Sep, 2014 1 commit
  25. 11 Sep, 2014 1 commit
    • Christopher Subich's avatar
      wave_reader, tracer file output · 2a7736a9
      Christopher Subich authored
      Previously for the wave_reader case file, "tracer.0" and
      "tracer_reader.m" files were being output regardless of whether the
      configuration specified that there was a passive tracer.  This arose
      because of a scoping error during tracer initalization; tracer.0
      contained a copy of rho.0.
      
      This has been fixed in this commit.  No numerical results were in error
      because of this bug.
      2a7736a9
  26. 10 Sep, 2014 1 commit
    • Michael Dunphy's avatar
      Embed the case file source in the SPINS binary. · 1d616e81
      Michael Dunphy authored
      The case file source code is embedded into the binary at compile time.
      At runtime, the source code is printed to spinscase.cpp. This eliminates
      the need to maintain a copy of the case file alongside the binary.
      1d616e81
  27. 17 Jun, 2014 1 commit
    • Christopher Subich's avatar
      Include sanity checking for quadrature weight init · c6a0c6d5
      Christopher Subich authored
      The quadrature weight "get" functions will now check whether the
      corresponding array has been initialized (length >= 1) before returning
      the array.  This will catch errors where quadrature weights are used
      before being initialized.
      c6a0c6d5
  28. 20 Mar, 2014 1 commit
    • Christopher Subich's avatar
      Fixed "automatic" grid calculation error for wave_reader · d2454796
      Christopher Subich authored
      The "automatic" (unmapped) grid calculation used by wave_reader
      (automatic_grid in BaseCase.cpp) was in error for grids of
      Chebyshev-type.  Instead of computing:
      
      z = -Lz/2 + Lz/2*cos(pi*(0:(N-1))/(N-1))
      
      it instead computed
      
      z = -Lz/2 + Lz/2*cos(pi*(0:(N-1))/N)
      
      (to use matlab-style notation).  This caused the written-out zgrid to be
      in error near the bottom (max z-index) boundary.  Files using a
      chebyshev-type unmapped x-coordinate would see the error as well in the
      x-grid (max x-index).
      
      This error affects only processing that uses the grid variables itself.
      Within wave_reader where this was used, the automatic-generated grids
      were used only for writing out; the input velocities and density were
      read directly in from their corresponding files.  Likewise, SPINS
      performed the "correct" differentiation internally without reference to
      the grid coordinate.
      d2454796
  29. 02 Oct, 2013 2 commits