fix hdf5 filewriting

and remove unnecessary functions in ResonanceInfo

sammy/src/endf/SammyResonanceInfo.h

@@ -110,55 +110,22 @@ namespace sammy{
         bool getIncludeInCalc() const{
             return includeInCalc;
         }
-        
-        /**
-         * Set the index of the covariance data for the resonance energy
-         * If -1, the covariance matrix does not contain information for this
-         * resonance energy
-         * 
-         * @param index index of the covariance data for the resonance energy
-         */
-        void setEnergyCovarianceIndex(int index){
-           if( (int)covarianceIndex.size() < 1) covarianceIndex.resize(1,-1);
-            covarianceIndex[0] = index;
-        }
-        
-        /**
-         * Get the index of the covariance data for the resonance energy
-         * If -1, the covariance matrix does not contain information for this
-         * resonance energy
-         * @return index of the covariance data for the resonance energy
-         */
-        int getEnergyCovarianceIndex() const{
-            if( (int)covarianceIndex.size() < 1) return -1;
-            return covarianceIndex[0];
-        }
-        
-        /**
-         * Set the index of the covariance data for the the indicated channel
-         * If -1, the covariance matrix does not contain information for this
-         * channel
-         * 
-         * @param index index of the covariance data for this channel
-         */
-        void setChannelCovarianceIndex(int nc, int index){
-           if( (int)covarianceIndex.size() < nc+2) covarianceIndex.resize(nc+2,-1);
-            covarianceIndex[nc+1] = index;
-        }
-        
-        /**
-         * Get the index of the covariance data for the the indicated channel
-         * If -1, the covariance matrix does not contain information for this
-         * resonance energy
-         * @return index of the covariance data for the resonance energy
-         */
-        int getChannelCovarianceIndex(int nc) const{
-            if( (int)covarianceIndex.size() < nc+2) return -1;
-            return covarianceIndex[nc+1];
-        }
-        
-        
        
+        /* TODO: These functions are doing double duty right now: they are used 
+           to set and get the "fit option" given in the input files, which is 
+           0,1, or 3 and they are used to get and set the covariance index. 
+
+           We should break this into two functions/variables so that we can 
+           expect when we will be getting/setting covariance indices or fitting
+           options.
+
+           Functions that need to change:
+           - setEnergyFitOption
+           - getEnergyFitOption
+           - setChannelFitOption
+           - getChannelFitOption
+        */
+
         /**
          * Set the fit  option  of the covariance data for the resonance energy.
          * 

sammy/src/endf/interface/cix/SammyResonanceInfo.cpp2f.xml

@@ -26,19 +26,6 @@
            <param name="incl"   type="bool"/> 
        </method>
        <method name="getIncludeInCalc" return_type="bool"/>
-        
-       <method name="setEnergyCovarianceIndex">
-           <param name="index"   type="int"  offset="-1"/> 
-       </method>
-       <method name="getEnergyCovarianceIndex" return_type="int"/>
-        
-       <method name="setChannelCovarianceIndex">
-           <param name="nc"   type="int"  offset="-1"/>
-           <param name="index"   type="int"  offset="-1"/> 
-       </method>
-       <method name="getChannelCovarianceIndex" return_type="int">
-           <param name="nc"   type="int"  offset="-1"/>  
-       </method>
        
         <method name="setEnergyFitOption">
            <param name="option"   type="int"/> 

sammy/src/endf/interface/cpp/SammyResonanceInfoInterface.cpp

@@ -2,7 +2,7 @@
 * This file has been dynamically generated by Class Interface Xml (CIX) 
 * DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 * If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Fri Feb 19 11:38:39 EST 2021
+* Date Generated: Tue Apr 05 15:00:41 EDT 2022
 * If any issues are experiences with this generated file that cannot be fixed
 * with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 */
@@ -54,26 +54,6 @@ bool SammyResonanceInfo_getIncludeInCalc(void * SammyResonanceInfo_ptr)
     return ((SammyResonanceInfo*)SammyResonanceInfo_ptr)->getIncludeInCalc();
 }
 
-void SammyResonanceInfo_setEnergyCovarianceIndex(void * SammyResonanceInfo_ptr,int * index)
-{
-    ((SammyResonanceInfo*)SammyResonanceInfo_ptr)->setEnergyCovarianceIndex(*index);
-}
-
-int SammyResonanceInfo_getEnergyCovarianceIndex(void * SammyResonanceInfo_ptr)
-{
-    return ((SammyResonanceInfo*)SammyResonanceInfo_ptr)->getEnergyCovarianceIndex();
-}
-
-void SammyResonanceInfo_setChannelCovarianceIndex(void * SammyResonanceInfo_ptr,int * nc,int * index)
-{
-    ((SammyResonanceInfo*)SammyResonanceInfo_ptr)->setChannelCovarianceIndex(*nc,*index);
-}
-
-int SammyResonanceInfo_getChannelCovarianceIndex(void * SammyResonanceInfo_ptr,int * nc)
-{
-    return ((SammyResonanceInfo*)SammyResonanceInfo_ptr)->getChannelCovarianceIndex(*nc);
-}
-
 void SammyResonanceInfo_setEnergyFitOption(void * SammyResonanceInfo_ptr,int * option)
 {
     ((SammyResonanceInfo*)SammyResonanceInfo_ptr)->setEnergyFitOption(*option);

sammy/src/endf/interface/cpp/SammyResonanceInfoInterface.h

@@ -2,7 +2,7 @@
 * This file has been dynamically generated by Class Interface Xml (CIX) 
 * DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 * If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Fri Feb 19 11:38:39 EST 2021
+* Date Generated: Tue Apr 05 15:00:41 EDT 2022
 * If any issues are experiences with this generated file that cannot be fixed
 * with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 */
@@ -22,10 +22,6 @@ void SammyResonanceInfo_setResonanceIndex(void * SammyResonanceInfo_ptr,int * in
 int SammyResonanceInfo_getResonanceIndex(void * SammyResonanceInfo_ptr);
 void SammyResonanceInfo_setIncludeInCalc(void * SammyResonanceInfo_ptr,bool * incl);
 bool SammyResonanceInfo_getIncludeInCalc(void * SammyResonanceInfo_ptr);
-void SammyResonanceInfo_setEnergyCovarianceIndex(void * SammyResonanceInfo_ptr,int * index);
-int SammyResonanceInfo_getEnergyCovarianceIndex(void * SammyResonanceInfo_ptr);
-void SammyResonanceInfo_setChannelCovarianceIndex(void * SammyResonanceInfo_ptr,int * nc,int * index);
-int SammyResonanceInfo_getChannelCovarianceIndex(void * SammyResonanceInfo_ptr,int * nc);
 void SammyResonanceInfo_setEnergyFitOption(void * SammyResonanceInfo_ptr,int * option);
 int SammyResonanceInfo_getEnergyFitOption(void * SammyResonanceInfo_ptr);
 void SammyResonanceInfo_setChannelFitOption(void * SammyResonanceInfo_ptr,int * nc,int * option);

sammy/src/endf/interface/fortran/SammyResonanceInfo_I.f90

@@ -2,7 +2,7 @@
 !! This file has been dynamically generated by Class Interface Xml (CIX) 
 !! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 !! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Fri Feb 19 11:38:39 EST 2021
+!! Date Generated: Tue Apr 05 15:00:41 EDT 2022
 !! If any issues are experiences with this generated file that cannot be fixed
 !! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 !!/
@@ -59,30 +59,6 @@ logical(C_BOOL) function f_SammyResonanceInfo_getIncludeInCalc(SammyResonanceInf
     implicit none
     type(C_PTR), value :: SammyResonanceInfo_ptr;
 end function
-subroutine f_SammyResonanceInfo_setEnergyCovarianceIndex(SammyResonanceInfo_ptr, index ) BIND(C,name="SammyResonanceInfo_setEnergyCovarianceIndex")
-    use,intrinsic :: ISO_C_BINDING
-    implicit none
-    type(C_PTR), value :: SammyResonanceInfo_ptr;
-    integer(C_INT) :: index;
-end subroutine
-integer(C_INT) function f_SammyResonanceInfo_getEnergyCovarianceIndex(SammyResonanceInfo_ptr ) BIND(C,name="SammyResonanceInfo_getEnergyCovarianceIndex")
-    use,intrinsic :: ISO_C_BINDING
-    implicit none
-    type(C_PTR), value :: SammyResonanceInfo_ptr;
-end function
-subroutine f_SammyResonanceInfo_setChannelCovarianceIndex(SammyResonanceInfo_ptr, nc,index ) BIND(C,name="SammyResonanceInfo_setChannelCovarianceIndex")
-    use,intrinsic :: ISO_C_BINDING
-    implicit none
-    type(C_PTR), value :: SammyResonanceInfo_ptr;
-    integer(C_INT) :: nc;
-    integer(C_INT) :: index;
-end subroutine
-integer(C_INT) function f_SammyResonanceInfo_getChannelCovarianceIndex(SammyResonanceInfo_ptr, nc ) BIND(C,name="SammyResonanceInfo_getChannelCovarianceIndex")
-    use,intrinsic :: ISO_C_BINDING
-    implicit none
-    type(C_PTR), value :: SammyResonanceInfo_ptr;
-    integer(C_INT) :: nc;
-end function
 subroutine f_SammyResonanceInfo_setEnergyFitOption(SammyResonanceInfo_ptr, option ) BIND(C,name="SammyResonanceInfo_setEnergyFitOption")
     use,intrinsic :: ISO_C_BINDING
     implicit none

sammy/src/endf/interface/fortran/SammyResonanceInfo_M.f90

@@ -2,7 +2,7 @@
 !! This file has been dynamically generated by Class Interface Xml (CIX) 
 !! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
 !! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Fri Feb 19 11:38:39 EST 2021
+!! Date Generated: Tue Apr 05 15:00:41 EDT 2022
 !! If any issues are experiences with this generated file that cannot be fixed
 !! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
 !!/
@@ -22,10 +22,6 @@ type SammyResonanceInfo
     procedure, pass(this) :: getResonanceIndex => SammyResonanceInfo_getResonanceIndex
     procedure, pass(this) :: setIncludeInCalc => SammyResonanceInfo_setIncludeInCalc
     procedure, pass(this) :: getIncludeInCalc => SammyResonanceInfo_getIncludeInCalc
-    procedure, pass(this) :: setEnergyCovarianceIndex => SammyResonanceInfo_setEnergyCovarianceIndex
-    procedure, pass(this) :: getEnergyCovarianceIndex => SammyResonanceInfo_getEnergyCovarianceIndex
-    procedure, pass(this) :: setChannelCovarianceIndex => SammyResonanceInfo_setChannelCovarianceIndex
-    procedure, pass(this) :: getChannelCovarianceIndex => SammyResonanceInfo_getChannelCovarianceIndex
     procedure, pass(this) :: setEnergyFitOption => SammyResonanceInfo_setEnergyFitOption
     procedure, pass(this) :: getEnergyFitOption => SammyResonanceInfo_getEnergyFitOption
     procedure, pass(this) :: setChannelFitOption => SammyResonanceInfo_setChannelFitOption
@@ -92,32 +88,6 @@ function SammyResonanceInfo_getIncludeInCalc(this) result(result2Return)
     logical(C_BOOL):: result2Return
     result2Return=f_SammyResonanceInfo_getIncludeInCalc(this%instance_ptr)
 end function
-subroutine SammyResonanceInfo_setEnergyCovarianceIndex(this, index)
-    implicit none
-    class(SammyResonanceInfo)::this
-    integer(C_INT)::index
-    call f_SammyResonanceInfo_setEnergyCovarianceIndex(this%instance_ptr, index-1)
-end subroutine
-function SammyResonanceInfo_getEnergyCovarianceIndex(this) result(result2Return)
-    implicit none
-    class(SammyResonanceInfo)::this
-    integer(C_INT):: result2Return
-    result2Return=f_SammyResonanceInfo_getEnergyCovarianceIndex(this%instance_ptr)  + 1
-end function
-subroutine SammyResonanceInfo_setChannelCovarianceIndex(this, nc, index)
-    implicit none
-    class(SammyResonanceInfo)::this
-    integer(C_INT)::nc
-    integer(C_INT)::index
-    call f_SammyResonanceInfo_setChannelCovarianceIndex(this%instance_ptr, nc-1,index-1)
-end subroutine
-function SammyResonanceInfo_getChannelCovarianceIndex(this, nc) result(result2Return)
-    implicit none
-    class(SammyResonanceInfo)::this
-    integer(C_INT)::nc
-    integer(C_INT):: result2Return
-    result2Return=f_SammyResonanceInfo_getChannelCovarianceIndex(this%instance_ptr, nc-1) + 1
-end function
 subroutine SammyResonanceInfo_setEnergyFitOption(this, option)
     implicit none
     class(SammyResonanceInfo)::this

sammy/src/io/Hdf5IO.cpp

@@ -103,7 +103,7 @@ namespace sammy {
 
     void Hdf5IO::writeRPCM(CovarianceData & cov,
                            SammyRMatrixParameters & resParData,
-                           std::string hdf5Filename,
+                           const std::string & hdf5Filename,
                            int rpcmIsReduced){
         if( rpcmIsReduced < 0 || rpcmIsReduced > 1 ){
             std::string errmsg("integer rpcmIsReduced must be 0 || 1");
@@ -113,23 +113,25 @@ namespace sammy {
         std::cout << "writing resonance covariance to HDF5 file..." << std::endl;
         std::cout << "--------------------------------------------" << std::endl;
 
+        // ---- set up the parameter with u-param to start ----------
         std::vector<double> param;
         int numParam = cov.getNumTotalParam();
         endf::ResonanceCovariance * parCov;
-        // ---- get parameter covariance and parameters out ---------
+        for( int i=0; i<numParam; ++i ) param.push_back(cov.getUParamValue(i));
+        // ---- get appropriate (p||u) covariance and parameters ----
         if( rpcmIsReduced ){
             parCov = cov.getUCovariance();
-            for( int i=0; i<numParam; ++i ) param.push_back(cov.getUParamValue(i));
         }
-        else{
+        else{ // get the physical pars
             CovToolBox ctb;
-            param = ctb.getPhysParam(resParData);
+            ctb.swapToPhysParam(resParData,param);
             parCov = cov.getCovariance();
         }
         // ---- check if length of params match covariance ----------
         int numRows = parCov->getNumRows();
-        if( numParam != numRows ) {
-            std::cerr << "numParam != numRows of covariance!" << std::endl;
+        if( numParam != numRows || numRows != param.size() ) {
+            std::string msg("covariance and parameters don't match in length!");
+            throw std::invalid_argument(msg);
         }
 
         double data[numRows][numRows];
@@ -148,13 +150,13 @@ namespace sammy {
         cdim[0] = numRows;
         cdim[1] = cdim[0];
         int crank = sizeof(cdim) / sizeof(hsize_t);
-        std::string covarname("/uCovar/");
+        std::string covarname("/covar/");
 
         // --- param space -----
         hsize_t pdim[1];
         pdim[0] = numParam;
         int prank = sizeof(pdim) / sizeof(hsize_t);
-        std::string paramname("/uParam/");
+        std::string paramname("/param/");
 
         H5::DataSpace cspace(crank,cdim);
         H5::DataSpace pspace(prank,pdim);
@@ -175,8 +177,8 @@ namespace sammy {
 
         try{
             H5::H5File file(hdf5Filename, H5F_ACC_RDONLY);
-            H5::DataSet h5covSet = file.openDataSet( "uCovar" );
-            H5::DataSet h5parSet = file.openDataSet( "uParam" );
+            H5::DataSet h5covSet = file.openDataSet( "covar" );
+            H5::DataSet h5parSet = file.openDataSet( "param" );
             H5::DataSpace h5covSpace = h5covSet.getSpace();
             H5::DataSpace h5parSpace = h5parSet.getSpace();
             // --- get dimensions of data to size local arrays

sammy/src/io/Hdf5IO.h

@@ -67,14 +67,14 @@ namespace sammy{
          ************************************************************/
         void writeRPCM(CovarianceData & cov,
                        SammyRMatrixParameters & resParData,
-                       std::string hdf5Filename,
+                       const std::string & hdf5Filename,
                        int rpcmIsReduced);
 
         /*************************************************************
          * Read the resonance parameter covariance matrix from an HDF5
          * file which was written by Hdf5IO::writeRPCM()
          * 
-         * @param cov class that we fill with param. covariance info
+         * @param cov class to be filled with param. covariance info
          * @param hdf5filename the string specifying the H5 file name
          ************************************************************/
         void readRPCM(CovarianceData & cov,std::string hdf5Filename);

sammy/src/salmon/CovToolBox.cpp

@@ -49,34 +49,42 @@ namespace sammy{
         }
     }
 
-    std::vector<double> CovToolBox::getPhysParam(SammyRMatrixParameters & resParData){
+    void CovToolBox::swapToPhysParam(SammyRMatrixParameters & resParData,
+                                  std::vector<double> & physParam){
 
         // --------------- grab all the *varied* pars --------------------------
+        int maxInd = physParam.size()-1;
         int nr = resParData.getNumResonances();
-        std::vector<double> physParam;
-        for(int iRes =0;iRes<nr;++iRes){ // loop over resonances
+        int covIndex = 0;
+        // loop over resonances
+        for(int iRes =0;iRes<nr;++iRes){
             sammy::SammyResonanceInfo *resInfo = resParData.getResonanceInfo(iRes);
             int iGrp = resInfo->getSpinGroupIndex();
             endf::RMatResonance *resonance = resParData.getResonance(resInfo);
             sammy::SammySpinGroupInfo *spinInfo = resParData.getSpinGroupInfo(iGrp);
-
             int parPerRes = spinInfo->getNumResPar();
-            for(int par=0;par<parPerRes;++par){ // loop over par's for each res
-                int covInd = 0; // >0 means varied or PUP'd
+
+            // loop over par's for each res
+            for(int par=0;par<parPerRes;++par){
+                int covInd = 0; // >=0 means varied or PUP'd
                 if(par==0){
-                    covInd = resInfo->getEnergyFitOption();
-                    if( covInd>0 && covInd!=3 ){
-                        physParam.push_back(resonance->getEres());
+                    covInd = resInfo->getEnergyFitOption()-1; //fortran indexed
+                    if( covInd>=0 && covInd<=maxInd ){
+                        physParam[covInd] = resonance->getEres();
                     }
                 }
                 else{
-                    covInd = resInfo->getChannelFitOption(par-1);
-                    if( covInd>0 && covInd!=3 ){
-                        physParam.push_back(resonance->getWidth(par-1));
+                    covInd = resInfo->getChannelFitOption(par-1)-1; //fortran indexed
+                    if( covInd>=0 && covInd<=maxInd ){
+                        physParam[covInd] = resonance->getWidth(par-1);
                     }
                 }
+                // throw if we have an indexing mistake
+                if( covInd > maxInd ) {
+                    std::string msg("covInd is incorrect for setting p-param!");
+                    throw std::runtime_error(msg);
+                }
             }
         }
-        return physParam;
     }
 }
\ No newline at end of file

sammy/src/salmon/CovToolBox.h

@@ -33,13 +33,15 @@ namespace sammy {
                                     endf::CrossCovList & cc);
 
         /***********************************************************************
-         * Find varied parameters fill a vector with them and return the 
-         * physical parameters in the rigid "SAMMY order" in a vector of doubles
+         * Find varied parameters and fill a vector of (typcially) u-params with
+         * the resonance physical parameters in the rigid "SAMMY order". Other 
+         * than resonance parameters, u-params should = p-params
          * 
          * @param resParData the full r-matrix parameters object
-         * @return physParam vector of varied physical parameters
+         * @param physParam vector of varied physical parameters
          **********************************************************************/
-        std::vector<double> getPhysParam(SammyRMatrixParameters & resParData);
+        void swapToPhysParam(SammyRMatrixParameters & resParData,
+                          std::vector<double> & physParam);
 
     private:
 

sammy/src/sam/sammy_api_m.f90

 module sammy_api_m
-    use, intrinsic :: ISO_C_BINDING
+    use, intrinsic :: iso_c_binding
     use, intrinsic :: iso_fortran_env ! input_unit
     use SammyFlowControl_M
     use over_common_m
@@ -182,11 +182,28 @@ contains
         call Prepare_Storage (Nnndat)
     end subroutine sammy_get_theory
 
+    subroutine sammy_calc_theory
+      if (urr_calc) then
+         call samacs_0
+         Where_To_Next = 'genfit'
+      else
+         CALL Samthe_0  ! 0K
+         call doppler_and_resolution_broadening() ! doppler and resolution broadening
+         Where_To_Next = 'genfit'
+         if (Kaverg.ne.0) Where_To_Next = 'samend'
+         IF (Kartgd.EQ.1)  Where_To_Next = 'samend'
+      end if
+    end subroutine
+
     subroutine sammy_advance_theory() BIND(C,name="sammy_advance_theory")
+        logical::moreData, moreIter
         if (Where_To_Next.eq.'samthe') return
 
         Where_To_Next = 'samfin'        
-        call sammy_do_segments('samthe')
+        ! call sammy_do_segments('samthe')
+        moreData = .false.
+        moreIter = .false.
+        call Samfin_0(moreData, moreIter)  ! convert RR parameters
     end subroutine sammy_advance_theory
 
     subroutine sammy_do_segments(endSeg)
@@ -287,19 +304,6 @@ contains
        Where_To_Next = 'samthe'
     end subroutine
 
-    subroutine sammy_calc_theory
-      if (urr_calc) then
-         call samacs_0
-         Where_To_Next = 'genfit'
-      else
-         CALL Samthe_0  ! 0K
-         call doppler_and_resolution_broadening() ! doppler and resolution broadening
-         Where_To_Next = 'genfit'
-         if (Kaverg.ne.0) Where_To_Next = 'samend'
-         IF (Kartgd.EQ.1)  Where_To_Next = 'samend'
-      end if
-    end subroutine
-
     subroutine sammy_next_segment()
         IF (Where_To_Next.EQ.Stop_Segment) THEN
             K_Stop = K_STOP + 1
@@ -322,6 +326,7 @@ contains
         ELSE
             WRITE (6,10000) Where_To_Next
 10000       FORMAT ('Where_To_Next is ', A6)
+            print *, "This is an unexpected route...gracefully failing."
             STOP
         END IF
     end subroutine sammy_next_segment