Merge pull request #653 from galacticusorg/featSigmaErrorMessage

Add a helpful error message when integration of σ(M) fails
galacticusorg · Jul 26, 2024 · 6020c20 · 6020c20
2 parents deaa21a + 8e58232
commit 6020c20
Showing 1 changed file with 22 additions and 4 deletions.
diff --git a/source/structure_formation.cosmological_mass_variance.filtered_power_spectrum.F90 b/source/structure_formation.cosmological_mass_variance.filtered_power_spectrum.F90
@@ -104,7 +104,8 @@ monotonic---for example in cases where $\sigma(M)$ becomes constant or close to
      ! Unique values in the variance table and their corresponding indices.
      type            (uniqueTable                            ), allocatable, dimension(:) :: rootVarianceUniqueTable
      logical                                                                              :: monotonicInterpolation                        , growthIsMassDependent_                               , &
-         &                                                                                   normalizationSigma8                   =.false., truncateAtParticleHorizon
+          &                                                                                  normalizationSigma8                   =.false., truncateAtParticleHorizon                            , &
+          &                                                                                  warnedNonIncreasing
    contains
      !![
      <methods>
@@ -370,6 +371,7 @@ function filteredPowerConstructorInternal(sigma8,amplitudeScalar,cosmologicalMas
        !!]
     end if
     self%initialized           =.false.
+    self%warnedNonIncreasing   =.false.
     self%growthIsMassDependent_=self%powerSpectrumPrimordialTransferred_%growthIsWavenumberDependent()
     self%fileName              =inputPath(pathTypeDataDynamic)                                                       // &
          &                      'largeScaleStructure/'                                                               // &
@@ -973,7 +975,10 @@ subroutine filteredPowerRetabulate(self,mass,time)
                    call Error_Report(message//{introspection:location})
                 else
                    message=message//char(10)//"         If problems occur consider not attempting to model structure below this mass scale."
-                   call Warn        (message                          )
+                   if (.not.self%warnedNonIncreasing) then
+                      call Warn     (message                          )
+                      self%warnedNonIncreasing=.true.
+                   end if
                 end if
              end if
           end do
@@ -994,11 +999,12 @@ double precision function rootVariance(time_,useTopHat)
       Compute the root-variance of mass in spheres enclosing the given {\normalfont \ttfamily mass} from the power spectrum.
       !!}
       use, intrinsic :: ISO_C_Binding           , only : c_size_t
-      use            :: Interface_GSL           , only : GSL_EBadTol      , GSL_ETol  , GSL_ERound, GSL_Success, &
+      use            :: Interface_GSL           , only : GSL_EBadTol      , GSL_ETol    , GSL_ERound , GSL_Success  , &
            &                                             GSL_EMaxIter
       use            :: Numerical_Constants_Math, only : Pi
       use            :: Numerical_Integration   , only : GSL_Integ_Gauss15, integrator
       use            :: Sorting                 , only : sort
+      use            :: Display                 , only : displayReset     , displayGreen, displayBlue, displayYellow
       implicit none
       double precision                         , intent(in   ) :: time_
       logical                                  , intent(in   ) :: useTopHat
@@ -1099,7 +1105,19 @@ double precision function rootVariance(time_,useTopHat)
          else if (status /= GSL_Success) then
             ! Integration failed for some other reason, report an error.
             computeLogarithmically=.false.
-            call Error_Report('integration over interval failed'//{introspection:location})
+            block
+              type     (varying_string) :: message
+              character(len=12        ) :: label
+              message='integration over interval failed'//char(10)//displayGreen()//'HELP:'//displayReset()//' try increasing the value of <'//displayBlue()//'tolerance'
+              if (useTopHat) then
+                 message=message//"TopHat"
+                 write (label,'(e12.6)') self%toleranceTopHat
+              else
+                 write (label,'(e12.6)') self%tolerance
+              end if
+              message=message//displayReset()//' '//displayYellow()//'value'//displayReset()//'='//displayGreen()//'"'//trim(adjustl(label))//'"'//displayReset()//'/> in the <'//displayBlue()//'cosmologicalMassVariance'//displayReset()//' '//displayYellow()//'value'//displayReset()//'='//displayGreen()//'"filteredPower"'//displayReset()//'> parameter'
+              call Error_Report(message//{introspection:location})
+            end block
          else if (integrandInterval <= 0.0d0 .and. wavenumberLower > 0.0d0) then
             ! Integration gave a zero result, and the lower limit is non-zero. This may be because the upper limit is large and
             ! the power is confined to small wavenumbers near the lower limit. This can happen, for example, if attempting to