text:
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION AccumulateLine(ID,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: AccumulateLine
        INTEGER(KIND=4)  :: AccumulateLineF
        AccumulateLine = AccumulateLineF(ID,LINE)
      END FUNCTION AccumulateLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION AddError(ID,ERROR_MSG)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: ERROR_MSG
        INTEGER(KIND=4)  :: AddError
        INTEGER(KIND=4)  :: AddErrorF
        AddError = AddErrorF(ID,ERROR_MSG)
      END FUNCTION AddError
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION AddWarning(ID,WARN_MSG)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: WARN_MSG
        INTEGER(KIND=4)  :: AddWarning
        INTEGER(KIND=4)  :: AddWarningF
        AddWarning = AddWarningF(ID,WARN_MSG)
      END FUNCTION AddWarning
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION ClearAccumulatedLines(ID)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: ClearAccumulatedLines
        INTEGER(KIND=4)  :: ClearAccumulatedLinesF
        ClearAccumulatedLines = ClearAccumulatedLinesF(ID)
      END FUNCTION ClearAccumulatedLines
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION CreateIPhreeqc()
        IMPLICIT NONE
        INTEGER(KIND=4)  :: CreateIPhreeqc
        INTEGER(KIND=4)  :: CreateIPhreeqcF
        CreateIPhreeqc = CreateIPhreeqcF()
      END FUNCTION CreateIPhreeqc
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION DestroyIPhreeqc(ID)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: DestroyIPhreeqc
        INTEGER(KIND=4)  :: DestroyIPhreeqcF
        DestroyIPhreeqc = DestroyIPhreeqcF(ID)
      END FUNCTION DestroyIPhreeqc
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetComponent(ID,N,COMP)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: COMP
        INTEGER(KIND=4)  :: GetComponent
        INTEGER(KIND=4)  :: GetComponentF
        GetComponent = GetComponentF(ID,N,COMP)
      END FUNCTION GetComponent
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetComponentCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetComponentCount
        INTEGER(KIND=4) :: GetComponentCountF
        GetComponentCount = GetComponentCountF(ID)
      END FUNCTION GetComponentCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetCurrentSelectedOutputUserNumber(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetCurrentSelectedOutputUserNumber
        INTEGER(KIND=4) :: GetCurrentSelectedOutputUserNumberF
        GetCurrentSelectedOutputUserNumber =
     &                     GetCurrentSelectedOutputUserNumberF(ID)
      END FUNCTION GetCurrentSelectedOutputUserNumber
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE GetDumpFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        CALL GetDumpFileNameF(ID,FNAME)
      END SUBROUTINE GetDumpFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetDumpFileOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetDumpFileOn
        INTEGER(KIND=4) :: GetDumpFileOnF
        IF (GetDumpFileOnF(ID).EQ.0) THEN
          GetDumpFileOn = .FALSE.
        ELSE
          GetDumpFileOn = .TRUE.
        ENDIF
      END FUNCTION GetDumpFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! GetDumpString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetDumpStringLine(ID,N,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: GetDumpStringLine
        INTEGER(KIND=4)  :: GetDumpStringLineF
        GetDumpStringLine = GetDumpStringLineF(ID,N,LINE)
      END FUNCTION GetDumpStringLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetDumpStringLineCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetDumpStringLineCount
        INTEGER(KIND=4) :: GetDumpStringLineCountF
        GetDumpStringLineCount = GetDumpStringLineCountF(ID)
      END FUNCTION GetDumpStringLineCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetDumpStringOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetDumpStringOn
        INTEGER(KIND=4) :: GetDumpStringOnF
        IF (GetDumpStringOnF(ID).EQ.0) THEN
          GetDumpStringOn = .FALSE.
        ELSE
          GetDumpStringOn = .TRUE.
        ENDIF
      END FUNCTION GetDumpStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE GetErrorFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        CALL GetErrorFileNameF(ID,FNAME)
      END SUBROUTINE GetErrorFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetErrorFileOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetErrorFileOn
        INTEGER(KIND=4) :: GetErrorFileOnF
        IF (GetErrorFileOnF(ID).EQ.0) THEN
          GetErrorFileOn = .FALSE.
        ELSE
          GetErrorFileOn = .TRUE.
        ENDIF
      END FUNCTION GetErrorFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetErrorOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetErrorOn
        INTEGER(KIND=4) :: GetErrorOnF
        IF (GetErrorOnF(ID).EQ.0) THEN
          GetErrorOn = .FALSE.
        ELSE
          GetErrorOn = .TRUE.
        ENDIF
      END FUNCTION GetErrorOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! GetErrorString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetErrorStringLine(ID,N,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: GetErrorStringLine
        INTEGER(KIND=4)  :: GetErrorStringLineF
        GetErrorStringLine = GetErrorStringLineF(ID,N,LINE)
      END FUNCTION GetErrorStringLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetErrorStringLineCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetErrorStringLineCount
        INTEGER(KIND=4) :: GetErrorStringLineCountF
        GetErrorStringLineCount = GetErrorStringLineCountF(ID)
      END FUNCTION GetErrorStringLineCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetErrorStringOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetErrorStringOn
        INTEGER(KIND=4) :: GetErrorStringOnF
        IF (GetErrorStringOnF(ID).EQ.0) THEN
          GetErrorStringOn = .FALSE.
        ELSE
          GetErrorStringOn = .TRUE.
        ENDIF
      END FUNCTION GetErrorStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE GetLogFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        CALL GetLogFileNameF(ID,FNAME)
      END SUBROUTINE GetLogFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetLogFileOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetLogFileOn
        INTEGER(KIND=4) :: GetLogFileOnF
        IF (GetLogFileOnF(ID).EQ.0) THEN
          GetLogFileOn = .FALSE.
        ELSE
          GetLogFileOn = .TRUE.
        ENDIF
      END FUNCTION GetLogFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! GetLogString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetLogStringLine(ID,N,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: GetLogStringLine
        INTEGER(KIND=4)  :: GetLogStringLineF
        GetLogStringLine = GetLogStringLineF(ID,N,LINE)
      END FUNCTION GetLogStringLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetLogStringLineCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetLogStringLineCount
        INTEGER(KIND=4) :: GetLogStringLineCountF
        GetLogStringLineCount = GetLogStringLineCountF(ID)
      END FUNCTION GetLogStringLineCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetNthSelectedOutputUserNumber(ID,N)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: N
        INTEGER(KIND=4) :: GetNthSelectedOutputUserNumber
        INTEGER(KIND=4) :: GetNthSelectedOutputUserNumberF
        GetNthSelectedOutputUserNumber =
     &                     GetNthSelectedOutputUserNumberF(ID,N)
      END FUNCTION GetNthSelectedOutputUserNumber
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetLogStringOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetLogStringOn
        INTEGER(KIND=4) :: GetLogStringOnF
        IF (GetLogStringOnF(ID).EQ.0) THEN
          GetLogStringOn = .FALSE.
        ELSE
          GetLogStringOn = .TRUE.
        ENDIF
      END FUNCTION GetLogStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE GetOutputFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        CALL GetOutputFileNameF(ID,FNAME)
      END SUBROUTINE GetOutputFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetOutputFileOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetOutputFileOn
        INTEGER(KIND=4) :: GetOutputFileOnF
        IF (GetOutputFileOnF(ID).EQ.0) THEN
          GetOutputFileOn = .FALSE.
        ELSE
          GetOutputFileOn = .TRUE.
        ENDIF
      END FUNCTION GetOutputFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! GetOutputString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetOutputStringLine(ID,N,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: GetOutputStringLine
        INTEGER(KIND=4)  :: GetOutputStringLineF
        GetOutputStringLine = GetOutputStringLineF(ID,N,LINE)
      END FUNCTION GetOutputStringLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetOutputStringLineCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetOutputStringLineCount
        INTEGER(KIND=4) :: GetOutputStringLineCountF
        GetOutputStringLineCount = GetOutputStringLineCountF(ID)
      END FUNCTION GetOutputStringLineCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetOutputStringOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetOutputStringOn
        INTEGER(KIND=4) :: GetOutputStringOnF
        IF (GetOutputStringOnF(ID).EQ.0) THEN
          GetOutputStringOn = .FALSE.
        ELSE
          GetOutputStringOn = .TRUE.
        ENDIF
      END FUNCTION GetOutputStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputColumnCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetSelectedOutputColumnCount
        INTEGER(KIND=4) :: GetSelectedOutputColumnCountF
        GetSelectedOutputColumnCount = GetSelectedOutputColumnCountF(ID)
      END FUNCTION GetSelectedOutputColumnCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetSelectedOutputCount
        INTEGER(KIND=4) :: GetSelectedOutputCountF
        GetSelectedOutputCount = GetSelectedOutputCountF(ID)
      END FUNCTION GetSelectedOutputCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE GetSelectedOutputFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        CALL GetSelectedOutputFileNameF(ID,FNAME)
      END SUBROUTINE GetSelectedOutputFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputFileOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetSelectedOutputFileOn
        INTEGER(KIND=4) :: GetSelectedOutputFileOnF
        IF (GetSelectedOutputFileOnF(ID).EQ.0) THEN
          GetSelectedOutputFileOn = .FALSE.
        ELSE
          GetSelectedOutputFileOn = .TRUE.
        ENDIF
      END FUNCTION GetSelectedOutputFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputRowCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetSelectedOutputRowCount
        INTEGER(KIND=4) :: GetSelectedOutputRowCountF
        GetSelectedOutputRowCount = GetSelectedOutputRowCountF(ID)
      END FUNCTION GetSelectedOutputRowCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputStringLine(ID,N,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: GetSelectedOutputStringLine
        INTEGER(KIND=4)  :: GetSelectedOutputStringLineF
        GetSelectedOutputStringLine =
     &                      GetSelectedOutputStringLineF(ID,N,LINE)
      END FUNCTION GetSelectedOutputStringLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputStringLineCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetSelectedOutputStringLineCount
        INTEGER(KIND=4) :: GetSelectedOutputStringLineCountF
        GetSelectedOutputStringLineCount =
     &                     GetSelectedOutputStringLineCountF(ID)
      END FUNCTION GetSelectedOutputStringLineCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputStringOn(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: GetSelectedOutputStringOn
        INTEGER(KIND=4) :: GetSelectedOutputStringOnF
        IF (GetSelectedOutputStringOnF(ID).EQ.0) THEN
          GetSelectedOutputStringOn = .FALSE.
        ELSE
          GetSelectedOutputStringOn = .TRUE.
        ENDIF
      END FUNCTION GetSelectedOutputStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetSelectedOutputValue(ID,ROW,COL,VTYPE,DVALUE,SVALUE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: ROW
        INTEGER(KIND=4)  :: COL
        INTEGER(KIND=4)  :: VTYPE
        REAL(KIND=8)     :: DVALUE
        CHARACTER(LEN=*) :: SVALUE
        INTEGER(KIND=4)  :: GetSelectedOutputValue
        INTEGER(KIND=4)  :: GetSelectedOutputValueF
        GetSelectedOutputValue = GetSelectedOutputValueF(ID,ROW,
     &                     COL,VTYPE,DVALUE,SVALUE)
      END FUNCTION GetSelectedOutputValue
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE GetVersionString(VERSION)
        IMPLICIT NONE
        CHARACTER(LEN=*) :: VERSION
        CALL GetVersionStringF(VERSION)
      END SUBROUTINE GetVersionString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! GetWarningString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetWarningStringLine(ID,N,LINE)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        CHARACTER(LEN=*) :: LINE
        INTEGER(KIND=4)  :: GetWarningStringLine
        INTEGER(KIND=4)  :: GetWarningStringLineF
        GetWarningStringLine = GetWarningStringLineF(ID,N,LINE)
      END FUNCTION GetWarningStringLine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION GetWarningStringLineCount(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTEGER(KIND=4) :: GetWarningStringLineCount
        INTEGER(KIND=4) :: GetWarningStringLineCountF
        GetWarningStringLineCount = GetWarningStringLineCountF(ID)
      END FUNCTION GetWarningStringLineCount
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION LoadDatabase(ID,FILENAME)
        IMPLICIT NONE
        INTEGER (KIND=4) :: ID
        CHARACTER(LEN=*) :: FILENAME
        INTEGER(KIND=4)  :: LoadDatabase
        INTEGER(KIND=4)  :: LoadDatabaseF
        LoadDatabase = LoadDatabaseF(ID,FILENAME)
      END FUNCTION LoadDatabase
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION LoadDatabaseString(ID,INPUT)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: INPUT
        INTEGER(KIND=4)  :: LoadDatabaseString
        INTEGER(KIND=4)  :: LoadDatabaseStringF
        LoadDatabaseString = LoadDatabaseStringF(ID,INPUT)
      END FUNCTION LoadDatabaseString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE OutputAccumulatedLines(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        CALL OutputAccumulatedLinesF(ID)
      END SUBROUTINE OutputAccumulatedLines
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE OutputErrorString(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        CALL OutputErrorStringF(ID)
      END SUBROUTINE OutputErrorString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      SUBROUTINE OutputWarningString(ID)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        CALL OutputWarningStringF(ID)
      END SUBROUTINE OutputWarningString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION RunAccumulated(ID)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: RunAccumulated
        INTEGER(KIND=4)  :: RunAccumulatedF
        RunAccumulated = RunAccumulatedF(ID)
      END FUNCTION RunAccumulated
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION RunFile(ID,FILENAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FILENAME
        INTEGER(KIND=4)  :: RunFile
        INTEGER(KIND=4)  :: RunFileF
        RunFile = RunFileF(ID,FILENAME)
      END FUNCTION RunFile
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION RunString(ID,INPUT)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: INPUT
        INTEGER(KIND=4)  :: RunString
        INTEGER(KIND=4)  :: RunStringF
        RunString = RunStringF(ID,INPUT)
      END FUNCTION RunString
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetBasicFortranCallback(ID,COOKIE)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        INTERFACE
          DOUBLE PRECISION FUNCTION cookie(x1, x2, str)
            DOUBLE PRECISION, INTENT(in) :: x1
            DOUBLE PRECISION, INTENT(in) :: x2
            CHARACTER(*), INTENT(in)     :: str
          END FUNCTION
        END INTERFACE
        INTEGER(KIND=4) :: SetBasicFortranCallback
        INTEGER(KIND=4) :: SetBasicFortranCallbackF
        SetBasicFortranCallback = SetBasicFortranCallbackF(ID,COOKIE)
      END FUNCTION SetBasicFortranCallback
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetCurrentSelectedOutputUserNumber(ID,N)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        INTEGER(KIND=4)  :: N
        INTEGER(KIND=4)  :: SetCurrentSelectedOutputUserNumber
        INTEGER(KIND=4)  :: SetCurrentSelectedOutputUserNumberF
        SetCurrentSelectedOutputUserNumber =
     &                      SetCurrentSelectedOutputUserNumberF(ID,N)
      END FUNCTION SetCurrentSelectedOutputUserNumber
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetDumpFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        INTEGER(KIND=4)  :: SetDumpFileName
        INTEGER(KIND=4)  :: SetDumpFileNameF
        SetDumpFileName = SetDumpFileNameF(ID,FNAME)
      END FUNCTION SetDumpFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetDumpFileOn(ID,DUMP_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: DUMP_ON
        INTEGER(KIND=4) :: SetDumpFileOn
        INTEGER(KIND=4) :: SetDumpFileOnF
        SetDumpFileOn = SetDumpFileOnF(ID,DUMP_ON)
      END FUNCTION SetDumpFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetDumpStringOn(ID,DUMP_STRING_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: DUMP_STRING_ON
        INTEGER(KIND=4) :: SetDumpStringOn
        INTEGER(KIND=4) :: SetDumpStringOnF
        SetDumpStringOn = SetDumpStringOnF(ID,DUMP_STRING_ON)
      END FUNCTION SetDumpStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetErrorFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        INTEGER(KIND=4)  :: SetErrorFileName
        INTEGER(KIND=4)  :: SetErrorFileNameF
        SetErrorFileName = SetErrorFileNameF(ID,FNAME)
      END FUNCTION SetErrorFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetErrorFileOn(ID,ERROR_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: ERROR_ON
        INTEGER(KIND=4) :: SetErrorFileOn
        INTEGER(KIND=4) :: SetErrorFileOnF
        SetErrorFileOn = SetErrorFileOnF(ID,ERROR_ON)
      END FUNCTION SetErrorFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetErrorOn(ID,ERROR_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: ERROR_ON
        INTEGER(KIND=4) :: SetErrorOn
        INTEGER(KIND=4) :: SetErrorOnF
        SetErrorOn = SetErrorOnF(ID,ERROR_ON)
      END FUNCTION SetErrorOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetErrorStringOn(ID,ERROR_STRING_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: ERROR_STRING_ON
        INTEGER(KIND=4) :: SetErrorStringOn
        INTEGER(KIND=4) :: SetErrorStringOnF
        SetErrorStringOn = SetErrorStringOnF(ID,ERROR_STRING_ON)
      END FUNCTION SetErrorStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetLogFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        INTEGER(KIND=4)  :: SetLogFileName
        INTEGER(KIND=4)  :: SetLogFileNameF
        SetLogFileName = SetLogFileNameF(ID,FNAME)
      END FUNCTION SetLogFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetLogFileOn(ID,LOG_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: LOG_ON
        INTEGER(KIND=4) :: SetLogFileOn
        INTEGER(KIND=4) :: SetLogFileOnF
        SetLogFileOn = SetLogFileOnF(ID,LOG_ON)
      END FUNCTION SetLogFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetLogStringOn(ID,LOG_STRING_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: LOG_STRING_ON
        INTEGER(KIND=4) :: SetLogStringOn
        INTEGER(KIND=4) :: SetLogStringOnF
        SetLogStringOn = SetLogStringOnF(ID,LOG_STRING_ON)
      END FUNCTION SetLogStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetOutputFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        INTEGER(KIND=4)  :: SetOutputFileName
        INTEGER(KIND=4)  :: SetOutputFileNameF
        SetOutputFileName = SetOutputFileNameF(ID,FNAME)
      END FUNCTION SetOutputFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetOutputFileOn(ID,OUTPUT_FILE_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: OUTPUT_FILE_ON
        INTEGER(KIND=4) :: SetOutputFileOn
        INTEGER(KIND=4) :: SetOutputFileOnF
        SetOutputFileOn = SetOutputFileOnF(ID,OUTPUT_FILE_ON)
      END FUNCTION SetOutputFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetOutputStringOn(ID,OUTPUT_STRING_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: OUTPUT_STRING_ON
        INTEGER(KIND=4) :: SetOutputStringOn
        INTEGER(KIND=4) :: SetOutputStringOnF
        SetOutputStringOn = SetOutputStringOnF(ID,OUTPUT_STRING_ON)
      END FUNCTION SetOutputStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetSelectedOutputFileName(ID,FNAME)
        IMPLICIT NONE
        INTEGER(KIND=4)  :: ID
        CHARACTER(LEN=*) :: FNAME
        INTEGER(KIND=4)  :: SetSelectedOutputFileName
        INTEGER(KIND=4)  :: SetSelectedOutputFileNameF
        SetSelectedOutputFileName = SetSelectedOutputFileNameF(ID,FNAME)
      END FUNCTION SetSelectedOutputFileName
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetSelectedOutputFileOn(ID,SELOUT_FILE_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: SELOUT_FILE_ON
        INTEGER(KIND=4) :: SetSelectedOutputFileOn
        INTEGER(KIND=4) :: SetSelectedOutputFileOnF
        SetSelectedOutputFileOn = SetSelectedOutputFileOnF(ID,
     &                     SELOUT_FILE_ON)
      END FUNCTION SetSelectedOutputFileOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      FUNCTION SetSelectedOutputStringOn(ID,SELOUT_STRING_ON)
        IMPLICIT NONE
        INTEGER(KIND=4) :: ID
        LOGICAL(KIND=4) :: SELOUT_STRING_ON
        INTEGER(KIND=4) :: SetSelectedOutputStringOn
        INTEGER(KIND=4) :: SetSelectedOutputStringOnF
        SetSelectedOutputStringOn = SetSelectedOutputStringOnF(ID,
     &                     SELOUT_STRING_ON)
      END FUNCTION SetSelectedOutputStringOn
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

ext: f

max_stars_repo_name: usgs-coupled/webmod

max_issues_repo_name: usgs-coupled/webmod

max_forks_repo_name: usgs-coupled/webmod


text:
MODULE m_statorbit

! ----------------------------------------------------------------------
! MODULE: m_statorbit.f03
! ----------------------------------------------------------------------
! Purpose:
!  Module for calling statorbit subroutine
! ----------------------------------------------------------------------
! Author :	Dr. Thomas Papanikolaou, Cooperative Research Centre for Spatial Information, Australia
! Created:	6 November 2017
! ----------------------------------------------------------------------

      IMPLICIT NONE
      !SAVE

Contains

SUBROUTINE statorbit (ds1, ds2, dorb_XYZ, dorb_RTN, dorb_Kepler, stat_XYZ, stat_RTN, stat_Kepler)

! ----------------------------------------------------------------------
! SUBROUTINE: statorbit.f03
! ----------------------------------------------------------------------
! Purpose:
!  Numerical comparison and statistics between two orbits in inertial frame (ICRF)
! ----------------------------------------------------------------------
! Input arguments:
! - orb_ref:		1st data series matrix (used as reference)
! - orb_det: 		2nd data series matrix
!
! Output arguments:
! - dorb_XYZ:		Orbit differences in the intertial frame (position and velocity vector differences)
!        			dorb_XYZ = [Time_MJD Time_Seconds00h r(X) r(Y) r(Z) v(X) v(Y) v(Z)]
! - dorb_RTN:		Orbit differences in the orbital frame (RTN: radial, along-track and cross-track differences)
!        			dorb_RTN = [Time_MJD Time_Seconds00h R(X) T(Y) N(Z) R(vX) T(vY) N(vZ)]
! - dorb_Kepler:	Orbit differences for the Keplerian elements
!        			dorb_Kepler = [Time_MJD Time_Seconds00h a e i Omega omega E]
!   				a: 		semi-major axis  (m)
!   				e: 		eccentricity
!   				i: 		inclination (degrees)
!   				Omega:	right ascension of the ascending node (degrees)
!   				omega:	argument of perigee (degrees)
!   				E: 		eccentric anomaly (degrees)
! - stat_XYZ:		Statistical quantities of dorb_XYZ (see Note 1)
! - stat_RTN:		Statistical quantities of dorb_RTN (see Note 1)
! - stat_Kepler:	Statistical quantities of dorb_Kepler (see Note 1)
! ----------------------------------------------------------------------
! Note 1
! The arrays of the statistical quantities of the orbital differences
! are formed as follows:
! stat_XYZ(1,:) : RMS value for each of the orbital elements e.g. X,Y,Z,Vx,Vy,Vz
! stat_XYZ(2,:) : Standard Deviation (sigma) for each of the orbital elements e.g. X,Y,Z,Vx,Vy,Vz
! stat_XYZ(3,:) : Mean value for each of the orbital elements e.g. X,Y,Z,Vx,Vy,Vz
! stat_XYZ(4,:) : MIN value for each of the orbital elements e.g. X,Y,Z,Vx,Vy,Vz
! stat_XYZ(5,:) : MAX value for each of the orbital elements e.g. X,Y,Z,Vx,Vy,Vz
! ----------------------------------------------------------------------
! Author :	Dr. Thomas Papanikolaou, Cooperative Research Centre for Spatial Information, Australia
! Created:	11 October 2017
! ----------------------------------------------------------------------

      USE mdl_precision
      USE mdl_num
      USE m_statist
      USE mdl_config
      IMPLICIT NONE

! ----------------------------------------------------------------------
! Dummy arguments declaration
! ----------------------------------------------------------------------
! IN
      REAL (KIND = prec_d), INTENT(IN), ALLOCATABLE, DIMENSION(:,:) :: ds1, ds2
! ----------------------------------------------------------------------
! OUT
      REAL (KIND = prec_d), INTENT(OUT), DIMENSION(:,:), ALLOCATABLE :: dorb_XYZ, dorb_RTN, dorb_Kepler
      ! REAL (KIND = prec_d), INTENT(OUT), DIMENSION(:), ALLOCATABLE :: RMSdsr, MINdsr, MAXdsr
      REAL (KIND = prec_d), INTENT(OUT), DIMENSION(5,6) :: stat_XYZ, stat_RTN, stat_Kepler
! ----------------------------------------------------------------------

! ----------------------------------------------------------------------
! Local variables declaration
! ----------------------------------------------------------------------
      REAL (KIND = prec_d), ALLOCATABLE, DIMENSION(:,:) :: dsr
      REAL (KIND = prec_d), ALLOCATABLE, DIMENSION(:) :: dx
      REAL (KIND = prec_d), DIMENSION(6) :: RMSdsr, Sigmadsr, MEANdsr, MINdsr, MAXdsr
      REAL (KIND = prec_d) :: RMSdx, Sigmadx, MEANdx, MINdx, MAXdx
! ----------------------------------------------------------------------
      REAL (KIND = prec_d) :: delta_t, dt_limit
      INTEGER (KIND = prec_int8) :: Nepochs, Nepochs2, Nelements, Nepochs_delta
      INTEGER (KIND = prec_int8) :: i, j, j1, k, looptest, l
      INTEGER (KIND = prec_int8) :: sz1, sz2, sz3, sz4
      INTEGER (KIND = prec_int2) :: AllocateStatus, DeAllocateStatus
! ----------------------------------------------------------------------
      REAL (KIND = prec_d), DIMENSION(3) :: r1i, v1i, r2i, v2i
      REAL (KIND = prec_d), DIMENSION(3) :: dr, dv, dr_RTN, dv_RTN
      REAL (KIND = prec_d) :: Rrtn(3,3)
      REAL (KIND = prec_d) :: GM
      REAL (KIND = prec_d) :: kepler1(6), kepler2(6), kepler_9(9)
      CHARACTER (LEN=100) :: mesg
      logical   allzero

      integer iy, im, id, jflag, ih, imin
      double precision    fd, fsec
! ----------------------------------------------------------------------

! Set the numerical differences digits limit for the common epochs
dt_limit = 1.D-08

! Gravitational constant
GM = 0.39860044150D+15

sz1 = size(ds1, DIM = 1)
sz2 = size(ds1, DIM = 2)
Nepochs = sz1

sz3 = size(ds2, DIM = 1)
sz4 = size(ds2, DIM = 2)
Nepochs2 = sz3
!print*,'statorbit: sz1, sz2 = ',sz1, sz2
!print*,'statorbit: sz3, sz4 = ',sz3, sz4
! ----------------------------------------------------------------------
! Test columns dimension
If (sz2 .NE. sz4) Then
write(mesg, *) "Input matrices dimension(DIM=2) do not agree. DIM = 2 (", &
        sz2, ",", sz4, ")"
call report ('WARNING', pgrm_name, 'statorbit', mesg, 'src/fortran/m_statorbit.f03', 1)
End If
! ----------------------------------------------------------------------

!debug: print out the integrated orbit values
if (.false.) then
    do j = 1, Nepochs2
        print *, "mjd = ", ds2(j,1) !mjd
        call iau_JD2CAL (2400000.5d0, ds2(j,1), Iy, Im, Id, Fd, JFLAG)
        Fd = fd * 86400.0d0
        ih = fd / 3600.0d0
        fd = fd - ih * 3600.0d0
        imin = fd / 60.0d0
        fsec = fd - imin * 60.0d0
        print *, Iy, Im, Id, ih, imin, fsec
        print *, "ICRF: ", ds2(j,3), ds2(j,4), ds2(j,5)
        ! TODO: convert to ITRF?
    end do
end if

! ----------------------------------------------------------------------
! Find the number of the common epochs
Nepochs_delta = 0
i = 0
j = 0
Do i = 1 , Nepochs
   ! Test the time argument:
   Do j = 1 , Nepochs2
      delta_t = ABS(ds2(j,1) - ds1(i,1))
      IF (delta_t < dt_limit) then
         ! Counter of the common epochs
         Nepochs_delta = Nepochs_delta + 1
         exit
      End if
   End Do
End Do
! ----------------------------------------------------------------------
!print *,"Nepochs_delta", Nepochs_delta
!print *,"Nepochs", Nepochs
!print *,"delta_t", delta_t

! Dynamic allocatable array
! Allocate the array of the numerical orbit comparison
ALLOCATE (dsr(Nepochs_delta,sz2), STAT = AllocateStatus)
if (AllocateStatus .ne. 0) then
        write (mesg, *) "failed to allocate dsr, dimensions = (", &
                Nepochs_delta, ",", sz2, ")"
call report ('FATAL', pgrm_name, 'statorbit', mesg, 'src/fortran/m_statorbit.f95', 1)
end if
ALLOCATE (dorb_XYZ(Nepochs_delta,sz2), STAT = AllocateStatus)
if (AllocateStatus .ne. 0) then
        write (mesg, *) "failed to allocate dorb_XYZ, dimensions = (", &
                Nepochs_delta, ",", sz2, ")"
call report ('FATAL', pgrm_name, 'statorbit', mesg, 'src/fortran/m_statorbit.f95', 1)
end if
ALLOCATE (dorb_RTN(Nepochs_delta,sz2), STAT = AllocateStatus)
if (AllocateStatus .ne. 0) then
        write (mesg, *) "failed to allocate dorb_RTN, dimensions = (", &
                Nepochs_delta, ",", sz2, ")"
call report ('FATAL', pgrm_name, 'statorbit', mesg, 'src/fortran/m_statorbit.f95', 1)
end if
ALLOCATE (dorb_Kepler(Nepochs_delta,sz2), STAT = AllocateStatus)
if (AllocateStatus .ne. 0) then
        write (mesg, *) "failed to allocate dorb_Kepler, dimensions = (", &
                Nepochs_delta, ",", sz2, ")"
call report ('FATAL', pgrm_name, 'statorbit', mesg, 'src/fortran/m_statorbit.f95', 1)
end if
!initialise
dsr = 0.d0
dorb_XYZ = 0.d0
dorb_RTN = 0.d0
dorb_kepler = 0.d0

! ----------------------------------------------------------------------
! Compute the numerical differences
i = 0
j = 0
k = 0
Do i = 1 , Nepochs
   ! Test the time argument:
   Do j = 1 , Nepochs2
      delta_t = ABS(ds2(j,1) - ds1(i,1))
      IF (delta_t < dt_limit) then
         k = k + 1
         if (.false.) then !debugging
         if (i==14 .and. j==74) then
             print *, "i = 14, j = 74, k = ", k
             print *, ds1(i, 3), ds1(i, 4), ds1(i, 5)
             print *, ds2(j, 3), ds2(j, 4), ds2(j, 5)
         end if
         end if !end debugging
	     ! State vector numerical differences at common epochs
		 dsr(k,1:2) = ds1(i,1:2)
         ! check ds1 and ds2 for allzeros - move on if so
         allzero = .true.
         do l = 3, sz2
             if (ds1(1, l) .ne. 0.0d0) then
                  allzero = .false.
             end if
         end do
         if (allzero) then
             k = k-1
             cycle
         end if
         allzero = .true.
         do l = 3, sz2
             if (ds2(j,l) .ne. 0.0d0) then
                     allzero = .false.
             endif
         end do
         if (allzero) then
             k = k-1
             cycle
         end if
         dsr(k,3:sz2) = ds2(j,3:sz2) - ds1(i,3:sz2)
IF(dsr(k,3)> 1.D+4.OR.dsr(k,4)> 1.D+4.OR.dsr(k,5)>1.D+4)THEN
PRINT*,'Something wrong with the numerical integration'
PRINT*,'Please adjust the integration step or maybe check your time system against the sp3 file being used'
print*,"Nepochs=",Nepochs, ", Nepochs2=", Nepochs2
print*, "previous values"
if (i > 1) then
        print*, ds1(i-1,3), ds1(i-1,4), ds1(i-1,5)
else
        print*, "no previous ds1 values"
endif
if (j > 1) then
        print*, ds2(j-1,3), ds2(j-1,4), ds2(j-1,5)
else
        print*, "no previous ds2 values"
end if
if (k>1) then
        print*, dsr(k-1,3), dsr(k-1,4), dsr(k-1,5)
else
        print*, "no previous difference values"
end if
PRINT*, "i=", i, ", j=", j, ", k=", k
print*, ds1(i,3), ds1(i,4), ds1(i,5)
print*, ds2(j,3), ds2(j,4), ds2(j,5)
print*, dsr(k,3), dsr(k,4), dsr(k,5)
print*, "next values"
if (i .ne. Nepochs) print*, ds1(i+1,3), ds1(i+1,4), ds1(i+1,5)
if (j .ne. Nepochs2) print*, ds2(j+1,3), ds2(j+1,4), ds2(j+1,5)
STOP
END IF
	     ! Orbital frame numerical differences at common epochs
		 r1i = ds1(i,3:5)
		 v1i = ds1(i,6:8)
		 dr = ds2(j,3:5) - ds1(i,3:5)
		 dv = ds2(j,6:8) - ds1(i,6:8)
		 CALL orb_frame(r1i, v1i, Rrtn)
		 CALL matrix_Rr(Rrtn, dr, dr_RTN)
		 CALL matrix_Rr(Rrtn, dv, dv_RTN)
		 dorb_RTN(k,1:2) = ds1(i,1:2)
         dorb_RTN(k,3:5) = dr_RTN
         dorb_RTN(k,6:8) = dv_RTN

	     ! Kepler elements numerical differences at common epochs
		 CALL kepler_z2k (r1i, v1i, GM, kepler_9)
		 kepler1(1:5) = kepler_9(1:5)
		 kepler1(6) = kepler_9(7)

		 r2i = ds2(j,3:5)
		 v2i = ds2(j,6:8)
		 CALL kepler_z2k (r2i, v2i, GM, kepler_9)
		 kepler2(1:5) = kepler_9(1:5)
		 kepler2(6) = kepler_9(7)

		 dorb_Kepler(k,1:2) = ds1(i,1:2)
         dorb_Kepler(k,3:8) = kepler2 - kepler1

	  End IF
   End Do
End Do
! ----------------------------------------------------------------------
!print *,"k", k

dorb_XYZ = dsr
sz1 = size(dsr, DIM = 1)

! ----------------------------------------------------------------------
! Statistics
! ----------------------------------------------------------------------

! Allocate arrays
ALLOCATE (dx(sz1), STAT = AllocateStatus)
if (AllocateStatus .ne. 0) then
        write (mesg, *) "failed to allocate dx, dimension = ", &
                Nepochs_delta
call report ('FATAL', pgrm_name, 'statorbit', mesg, 'src/fortran/m_statorbit.f95', 1)
end if
dx = 0.d0

! Inertial frame
dsr = dorb_XYZ
i = 0
j = 0
Do i = 3 , sz2
   dx(:) = dsr(:,i)
   Call statist(dx, RMSdx, Sigmadx, MEANdx, MINdx, MAXdx)
   j = j + 1
   RMSdsr(j) = RMSdx
   MINdsr(j) = MINdx
   MAXdsr(j) = MAXdx
   Sigmadsr(j) = Sigmadx
   MEANdsr(j) = MEANdx
End Do
stat_XYZ(1,:) = RMSdsr
stat_XYZ(2,:) = Sigmadsr
stat_XYZ(3,:) = MEANdsr
stat_XYZ(4,:) = MINdsr
stat_XYZ(5,:) = MAXdsr

! Orbital frame
dsr = dorb_RTN
i = 0
j = 0
Do i = 3 , sz2
   dx(:) = dsr(:,i)
   Call statist(dx, RMSdx, Sigmadx, MEANdx, MINdx, MAXdx)
   j = j + 1
   RMSdsr(j) = RMSdx
   MINdsr(j) = MINdx
   MAXdsr(j) = MAXdx
   Sigmadsr(j) = Sigmadx
   MEANdsr(j) = MEANdx
End Do
stat_RTN(1,:) = RMSdsr
stat_RTN(2,:) = Sigmadsr
stat_RTN(3,:) = MEANdsr
stat_RTN(4,:) = MINdsr
stat_RTN(5,:) = MAXdsr

! Kepler frame
dsr = dorb_Kepler
i = 0
j = 0
Do i = 3 , sz2
   dx(:) = dsr(:,i)
   Call statist(dx, RMSdx, Sigmadx, MEANdx, MINdx, MAXdx)
   j = j + 1
   RMSdsr(j) = RMSdx
   MINdsr(j) = MINdx
   MAXdsr(j) = MAXdx
   Sigmadsr(j) = Sigmadx
   MEANdsr(j) = MEANdx
End Do
stat_Kepler(1,:) = RMSdsr
stat_Kepler(2,:) = Sigmadsr
stat_Kepler(3,:) = MEANdsr
stat_Kepler(4,:) = MINdsr
stat_Kepler(5,:) = MAXdsr
! ----------------------------------------------------------------------

 100 END SUBROUTINE

End Module

max_stars_repo_name: HiTMonitor/ginan

max_issues_repo_name: HiTMonitor/ginan

max_forks_repo_name: HiTMonitor/ginan


text:
The Human Powered Housing Project is a bike/RV contraption built and lived in by Brian Campbell. Its about fifteen feet long and is absolutely one of a kind. Brian rides around in style under a shade canopy while protected by a motorcycle like windshield. His house forms the back of the rig. It has windows, a slot for donations, a comfortable looking bed, shelves and everything else a normal person would want. When parked, the project sits on an incredibly sturdy kick stand because this is probably the heaviest bike that youve ever seen. Its made mostly out of materials that you can find at a hardware store. It looks like there are a lot of aluminum struts and foil tape. The back wheel is much larger than your average bike wheel but the front, while still homemade, looks more normal. The gear system is pretty beefy as well. This is at least the second Human Powered Housing Project to exist. Look below for a picture of an earlier prototype. Anyone know what this thing weighs?

Allegedly, he rents out his room for $10/night. Thats probably less than what youre paying for rent. Rooms? How many does it have? Has anyone ever done so?

http://bikeportland.org/2009/08/24/avisitfrommotorhomebikebuilderbriancampbell22672 Brian is now living in the Portland Oregon area.  He hasnt lived in Davis for a few years. My source, who was formerly a homeless street musician, now a bike builder, said that about five years ago he verified Brians riding when they both rode up Siskiyou Mountain Pass.

http://bullsheet.wordpress.com/2009/01/13/sleepingbiketrailers/ Sleeping bike trailers have been around for many years. There are others who live on their bike, in a http://www.bikeroute.com/michaelwashington/ different way.

On a related subject, Davis is world famous for its Bicycles humanpowered vehicles. There are many bikes that are WhymCycles unique to Davis.

The http://microship.com/resources/winnebikobehemoth.html Winnebiko of hightech nomad Steve Roberts was an earlier version of an allinone bicyclecentered environment.

As seen in the http://www.californiaaggie.com/media/storage/paper981/news/2003/11/14/FrontPage/Pedal.To.The.Metal1315811.shtml?norewrite200604171650&sourcedomainwww.californiaaggie.com Aggie.

 Questions
Brian claims to be able to sustain 7080 mph on the vehicle.  He has told me how he does it, but I have not seen him go anywhere near that.  Anyone know for a fact whether it is possible, and if so, whether he actually does it? Users/ChristopherMcKenzie

  7080 mph?  Perhaps this may be the vehicles terminal velocity after being ridden off the edge of a cliff  ;)  More seriously, the top speed in a downhillcoast situation depends crucially on the vehicles wiki:wikipedia:drag coefficient (air friction).  And the more drag there is, the more weight becomes a factor (heavierfaster).  The Human Powered Housing Project is neither sleek nor heavy  I would not bet on it in a wiki:wikipedia:Soap Box Derby, but I would wager that in a windtunnel situation at 80mph, all that the foil tape you see would lose its grip on those foam panels!  Users/ErikAnderson E
  That is when his front cowcatcher gets caught on the rear bumper of an 18wheeler going down I5.
    Ive seen a childs car which goes that fast. The speedometer was painted on though...

Doesnt the rear cabin of this bike look a lot shorter now than when this bike originally appeared in Davis? I think it may have been recently truncated? Was there some sort of damage to the rear which was removed? It doesnt seem like a good nights sleep can be had in it now.Users/RickEle

 Photos

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


!>
!> @file
!>
!>  PLASMA is a software package provided by:
!>  University of Tennessee, US,
!>  University of Manchester, UK.
!>
!> @precisions normal z -> s d c
!>
program test_zgeqrf_omp
    use, intrinsic :: iso_fortran_env

    use plasma
    implicit none

    integer, parameter :: sp = REAL32
    integer, parameter :: dp = REAL64

    ! set working precision, this value is rewritten for different precisions
    integer, parameter :: wp = dp

    ! matrices
    complex(wp), allocatable :: A(:,:)
    complex(wp), allocatable :: Aref(:,:)

    ! matrix descriptor for storing coefficients of the Householder reflectors
    ! generalization of TAU array known from LAPACK
    type(plasma_desc_t) :: descT

    ! descriptor for A
    type(plasma_desc_t) :: descA

    ! matrices for checking correctness
    complex(wp), allocatable :: Q(:,:)
    complex(wp), allocatable :: Id(:,:)
    complex(wp), allocatable :: R(:,:)

    ! dimensions
    integer :: m
    integer :: n
    integer :: lda
    integer :: nb
    integer :: ib
    logical :: tree_householder
    logical :: use_tuning

    real(wp) :: tol

    ! auxiliary variables
    integer :: ldq
    integer :: info
    integer :: minmn

    complex(wp), parameter :: zone  =  1.0
    complex(wp), parameter :: zmone = -1.0
    complex(wp), parameter :: zzero =  0.0
    integer :: seed(4) = [0, 0, 0, 1]

    real(wp) :: normA
    real(wp) :: ortho
    real(wp) :: error

    real(wp), allocatable :: work(:)
    real(wp) :: zlange, zlanhe, dlamch

    integer :: householder_mode
    integer :: householder_mode_check
    integer :: tuning_mode

    type(plasma_context_t), pointer :: plasma_context
    integer(kind=c_size_t) :: lworkspace
    type(plasma_workspace_t) :: workspace
    type(plasma_sequence_t) :: plasma_sequence
    type(plasma_request_t) :: plasma_request

    integer :: plasma_status

    ! Set parameters.
    m    = 2000 ! number of rows
    n    = 500  ! number of columns

    nb   = 256  ! tile size (square tiles)
    ib   = 64   ! inner blocking size within a tile

    ! Should tree reduction be used? This is good for tall and skinny matrices.
    tree_householder = .true.

    ! Should tuning be used?
    ! If it is on, the tests below reset the nb and ib values and then fail.
    use_tuning = .false.

    ! tolerance for the correctness check
    tol = 50.0 * dlamch('E')

    ! Initialize arrays.
    lda  = m    ! leading dimension of A
    allocate(A(lda, n), Aref(lda, n))
    call zlarnv(1, seed, lda*n, A)
    Aref = A

    write(*,'(a, i6, a, i6)') " Testing asynchronous QR factorization of a matrix ", m, " x ", n
    write(*,'(a, i6)')        "  tile size:        ", nb
    write(*,'(a, i6)')        "  inner block size: ", ib
    if (tree_householder) then
        write(*,'(a, i6)')    "  Using tree-based Householder reduction. "
    else
        write(*,'(a, i6)')    "  Using flat-tree Householder reduction. "
    end if
    if (use_tuning) then
        write(*,'(a, i6)')    "  Tuning enabled. "
    else
        write(*,'(a, i6)')    "  Tuning disabled. "
    end if

    !==============================================
    ! Initialize PLASMA.
    !==============================================
    call plasma_init(info)
    call check_error('plasma_init()', info)

    !==============================================
    ! Set PLASMA parameters.
    !==============================================
    ! set tile size
    call plasma_set(PlasmaNb, nb, info)
    call check_error('plasma_set()', info)

    ! set inner blocking size
    call plasma_set(PlasmaIb, ib, info)
    call check_error('plasma_set()', info)

    ! set Householder mode
    if (tree_householder) then
        householder_mode = PlasmaTreeHouseholder
    else
        householder_mode = PlasmaFlatHouseholder
    end if
    call plasma_set(PlasmaHouseholderMode, householder_mode, info)
    call check_error('plasma_set()', info)

    call plasma_get(PlasmaHouseholderMode, householder_mode_check, info)
    if (householder_mode_check /= householder_mode) then
        write(*,'(a, i2, a, i2)') &
            "Something is wrong - Householder mode not set correctly:", &
            householder_mode_check, " not the same as ", householder_mode
    end if

    ! set tuning
    if (use_tuning) then
        tuning_mode = PlasmaEnabled
    else
        tuning_mode = PlasmaDisabled
    end if
    call plasma_set(PlasmaTuning, tuning_mode, info)
    call check_error('plasma_set()', info)

    ! get PLASMA context
    call plasma_context_self(plasma_context)
    call check_error('plasma_context_self()', info)

    ! set tuning
    if (tuning_mode == PlasmaEnabled) then
        call plasma_tune_geqrf(plasma_context, PlasmaComplexDouble, m, n)
        call check_error('plasma_tune_geqrf()', info)
    end if

    ! create descriptor for A
    call plasma_desc_general_create(PlasmaComplexDouble, nb, nb, &
                                    m, n, 0, 0, m, n, descA, info)
    call check_error('plasma_desc_general_create()', info)

    ! create descriptor for T
    call plasma_descT_create(descA, ib, householder_mode, descT, info)
    call check_error('plasma_descT_create()', info)

    ! Allocate workspace.
    lworkspace = nb + ib*nb;  ! geqrt: tau + work
    call plasma_workspace_create(workspace, lworkspace, PlasmaComplexDouble, info)
    call check_error('plasma_workspace_create()', info)

    ! Initialize sequence.
    call plasma_sequence_init(plasma_sequence, info)
    call check_error('plasma_sequence_init()', info)

    ! Initialize request.
    call plasma_request_init(plasma_request, info)
    call check_error('plasma_request_init()', info)

    !==============================================
    ! Call QR factorization for A.
    !==============================================
    ! asynchronous block
    !$omp parallel
    !$omp master
        ! Translate to tile layout.
        call plasma_omp_zge2desc(A, lda, descA, plasma_sequence, &
                                 plasma_request)

        ! Call the tile async function.
        call plasma_omp_zgeqrf(descA, descT, workspace, plasma_sequence, &
                               plasma_request)

        ! Translate back to LAPACK layout.
        call plasma_omp_zdesc2ge(descA, A, lda, plasma_sequence, &
                                 plasma_request)
    !$omp end master
    !$omp end parallel
    ! implicit synchronization

    plasma_status = plasma_sequence%status
    print *, "PLASMA status after the OMP block: ", plasma_status

    call plasma_workspace_destroy(workspace, info)
    call check_error('plasma_workspace_destroy()', info)

    ! Free matrix A in tile layout.
    call plasma_desc_destroy(descA, info)
    call check_error('plasma_desc_destroy()', info)

    ! Test results by checking orthogonality of Q and precision of Q*R.
    ! Allocate space for Q.
    minmn = min(m, n)
    ldq = m
    allocate(Q(ldq, minmn))

    !==================================================================
    ! Build Q.
    !==================================================================
    call plasma_zungqr(m, minmn, minmn, A, lda, descT, Q, ldq, info)
    call check_error('plasma_zungqr()', info)

    ! Build the identity matrix.
    allocate(Id(minmn, minmn))
    call zlaset('g', minmn, minmn, zzero, zone, Id, minmn)

    ! Perform Id - Q^H * Q
    call zherk('Upper', 'Conjugate transpose', minmn, m, &
               zmone, Q, ldq, zone, Id, minmn)

    ! |Id - Q^H * Q|_oo
    allocate(work(m))
    ortho = zlanhe('I', 'u', minmn, Id, minmn, work)
    print *, 'ortho', ortho

    ! normalize the result
    ! |Id - Q^H * Q|_oo / n
    ortho = ortho / minmn

    deallocate(Q)
    deallocate(Id)

    ! Check the accuracy of A - Q * R
    ! Extract the R factor.
    allocate(R(m,n))
    R = 0.0
    call zlacpy('u', m, n, A, lda, R, m)

    !==================================================================
    ! Compute Q * R.
    !==================================================================
    call plasma_zunmqr(PlasmaLeft, PlasmaNoTrans, m, n, minmn, A, lda, descT, &
                       R, m, info)
    call check_error('plasma_zunmqr()', info)

    ! Compute the difference.
    ! R = A - Q*R
    R = Aref(1:m,:) - R

    ! |A|_oo
    normA = zlange('I', m, n, Aref, lda, work)

    ! |A - Q*R|_oo
    error = zlange('I', m, n, R, m, work)

    ! normalize the result
    ! |A-QR|_oo / (|A|_oo * n)
    error = error / (normA * n)

    ! print output of the check
    write(*,'(a, e13.6)') &
        "  Orthogonality of Q, measured as |Id - Q^H * Q|_oo / n: ", ortho
    write(*,'(a, e13.6)') &
        "  Error of factorization, measured as |A-Q*R|_oo / (|A|_oo * n): ", error
    if (error < tol .and. ortho < tol) then
        write(*,'(a)') "  The result is CORRECT."
    else
        write(*,'(a)') "  The result is WRONG!"
    end if
    write(*,'(a)') ""

    deallocate(work)
    deallocate(R)

    !==============================================
    ! Deallocate matrix in descriptor.
    !==============================================
    call plasma_desc_destroy(descT, info)
    call check_error('plasma_desc_destroy()', info)

    !==============================================
    ! Finalize PLASMA.
    !==============================================
    call plasma_finalize(info)
    call check_error('plasma_finalize()', info)

    deallocate(A)
    deallocate(Aref)

contains

    subroutine check_error(routine_name, code)
        use plasma
        implicit none
        character(len=*) :: routine_name
        integer :: code

        if (code /= PlasmaSuccess) then
            write(*,'(a, i2, a, a)') "Error ", code, "in a call to ", &
                trim(routine_name)
        end if
    end subroutine

end program test_zgeqrf_omp


! Copyright (c) 2019, NVIDIA CORPORATION.  All rights reserved.
!
! Licensed under the Apache License, Version 2.0 (the "License");
! you may not use this file except in compliance with the License.
! You may obtain a copy of the License at
!
!     http://www.apache.org/licenses/LICENSE-2.0
!
! Unless required by applicable law or agreed to in writing, software
! distributed under the License is distributed on an "AS IS" BASIS,
! WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
! See the License for the specific language governing permissions and
! limitations under the License.

module m
  type :: t
    procedure(a), pointer, pass :: c
    procedure(a), pointer, pass(x) :: d
  contains
    procedure, pass(y) :: a, b
  end type
contains
  subroutine a(x, y)
    class(t) :: x, y
  end
  subroutine b(y)
    class(t) :: y
  end
end module

!Expect: m.mod
!module m
!  type::t
!    procedure(a),pass,pointer::c
!    procedure(a),pass(x),pointer::d
!  contains
!    procedure,pass(y)::a
!    procedure,pass(y)::b
!  end type
!contains
!  subroutine a(x,y)
!    class(t)::x
!    class(t)::y
!  end
!  subroutine b(y)
!    class(t)::y
!  end
!end


text:
program run_dStar
	use iso_fortran_env, only : output_unit, error_unit
    use NScool_def
    use NScool_lib
    use constants_def, only : boltzmann
    use alt_micro
    use argparse

    character(len=*), parameter :: default_dStar_dir = '../../dStar'
    character(len=*), parameter :: default_inlist_file = 'inlist'
    character(len=64) :: my_dStar_dir
    character(len=64) :: inlist
    real(dp) :: eV_to_MK
    type(NScool_info), pointer :: s=>null()
    integer :: ierr, NScool_id, i
    real(dp), dimension(8) :: pred_Teff, obs_Teff, obs_Teff_sig
    real(dp) :: chi2

    ierr = 0
    call command_arg_set( &
        & 'dStar_directory',"sets the main dStar root directory",ierr, &
        & flag='D',takes_parameter=.TRUE.)
    call check_okay('set command argument dStar_directory',ierr)

    call command_arg_set( &
        & 'inlist_file','sets the namelist parameter file',ierr, &
        & flag='I',takes_parameter=.TRUE.)
    call check_okay('set command argument inlist file',ierr)

    call parse_arguments(ierr)
    call check_okay('parse_arguments',ierr)

    my_dStar_dir = trim(get_command_arg('dStar_directory'))
    if (len_trim(my_dStar_dir)==0) my_dStar_dir = default_dStar_dir
    inlist = trim(get_command_arg('inlist_file'))
    if (len_trim(inlist)==0) inlist = default_inlist_file

    call NScool_init(my_dStar_dir, ierr)
    call check_okay('NScool_init',ierr)

    NScool_id = alloc_NScool(ierr)
    call check_okay('NScool_id',ierr)

    call NScool_setup(NScool_id,inlist,ierr)
    call check_okay('NScool_setup',ierr)

    call get_NScool_info_ptr(NScool_id,s,ierr)
    call check_okay('get_NScool_info_ptr',ierr)

    s% other_set_Qimp => alt_Qimp
    s% other_sf_get_results => alt_sf

    call NScool_create_model(NScool_id,ierr)
    call check_okay('NScool_create_model',ierr)

    call NScool_evolve_model(NScool_id,ierr)
    call check_okay('NScool_evolve_model',ierr)

    ! we don't want to compare the effective temp. at t = 0, the end of the
    ! outburst
    pred_Teff = s% Teff_monitor(2:)/1.0e6
    eV_to_MK = 1.602176565e-12_dp/boltzmann/1.0e6
    ! observed effective temperatures (eV) and uncertainties
    obs_Teff = [103.2,88.9,75.5,73.3,71.0,66.0,70.3,63.1] * eV_to_MK
    obs_Teff_sig = [1.7,1.3,2.2,2.3,1.8,4.5,2.1,2.1] * eV_to_MK
    chi2 = sum((pred_Teff-obs_Teff)**2 / obs_Teff_sig**2 )

    write(output_unit,*)
    write(output_unit,'(a7,a6,tr3,a12)') 'time','Teff','obs. range'
    write(output_unit,'(a7,a6,tr3,a12)') '[d]','[MK]','[MK]'
    write(output_unit,'(13("-"),tr3,12("-"))')
    do i = 1, 8
        write(output_unit,'(f7.1,f6.3,tr3,f6.3,f6.3)')  &
        & s% t_monitor(i+1), &
        & pred_Teff(i),(obs_Teff(i)-obs_Teff_sig(i)), &
        & (obs_Teff(i)+obs_Teff_sig(i))
    end do
    write(output_unit,*)
    write(output_unit,'(a,f6.2)') 'chi2 = ',chi2

    call NScool_shutdown

contains
	subroutine check_okay(msg,ierr)
		character(len=*), intent(in) :: msg
		integer, intent(inout) :: ierr
		if (ierr /= 0) then
			write (error_unit,*) trim(msg)//': ierr = ',ierr
			if (ierr < 0) stop
		end if
	end subroutine check_okay
end program run_dStar

max_stars_repo_name: LBJ-Wade/dStar

max_issues_repo_name: LBJ-Wade/dStar

max_forks_repo_name: LBJ-Wade/dStar


text:
!-----------------------------------------------------------------------
!    Module:       prb_cntl_module
!    Author:       S. W. Bruenn
!    Date:         8/16/02
!-----------------------------------------------------------------------

MODULE prb_cntl_module

USE kind_module

!-----------------------------------------------------------------------
!
!             \\\\\ ZERO TRANSVERSE VELOCITY OPTION /////
!
!-----------------------------------------------------------------------
!  v_trans_0      : zero transverse velocity above shock switch
!
!     v_trans_0 = ye : transverse velocities 0 above shock
!     v_trans_0 = no : transverse velocities 0 above shock computed
!-----------------------------------------------------------------------

CHARACTER(LEN=2)                           :: v_trans_0

!-----------------------------------------------------------------------
!  Hydrodynamic Controls
!-----------------------------------------------------------------------
!  ihydro   : Hydrodynamics switch.
!
!     ihydro = 0   : hydrodynamics bypassed.
!     ihydro = 1   : hydrodynamics included.
!
!  irelhy   : Newtonian - GR hydrodynamics toggle.
!
!     irelhy = 0   : Newtonian hydrodynamics.
!     irelhy = 1   : PN hydrodynamics.
!     irelhy = 2   : general relativistic hydrodynamics.
!
!  ilapsehy : Newtonian - PN - GR lapse toggle.
!
!     ilapsehy = 0 : lapse = 1 in the hydrodynamics and transport.
!     ilapsehy = 1 : lapse computed to PN approximation in the /; b4..040
!      hydrodynamics if irelhy = 1, and in the transport if ireltrns = 1.
!     ilapsehy = 2 : lapse computed with full GR in the
!      hydrodynamics if irelhy = 1, and in the transport if ireltrns = 1.
!
!-----------------------------------------------------------------------

INTEGER                                    :: ihydro
INTEGER                                    :: irelhy
INTEGER                                    :: ilapsehy

!-----------------------------------------------------------------------
!  Hydrodynamic Mode Control
!-----------------------------------------------------------------------
!
!     e_mode = 'e  ' : evolve total energy (internal + gravitational + kinetic)
!     e_mode = 'e-g' : evolve partial energy (internal + kinetic)
!-----------------------------------------------------------------------

CHARACTER(LEN=3)                           :: e_mode

!-----------------------------------------------------------------------
!  Gravitation Options
!-----------------------------------------------------------------------
!  i_grav     : gravitation key
!
!     i_grav = 1 : spherical symmetric gravity
!     i_grav = 2 : nonspherical components computed from by a Newtonian
!                   Poisson solver
!
!-----------------------------------------------------------------------

INTEGER                                    :: i_grav

!-----------------------------------------------------------------------
!
!                     \\\\\ MOVING GRID OPTION /////
!
!-----------------------------------------------------------------------
!  m_grid_1 = la : Radial grid above the inner-outer grid boundary
!                   is Lagrangian if lagr and eul = no. (Control)
!  m_grid_1 = no : Radial grid above the inner-outer grid boundary
!                   is Eulerian if lagr and eul = no. (Control)
!
!  m_grid_2       This option is not currently functional
!
!-----------------------------------------------------------------------

CHARACTER(LEN=2)                           :: m_grid_1
CHARACTER(LEN=2)                           :: m_grid_2

!-----------------------------------------------------------------------
!  Transport controls
!-----------------------------------------------------------------------
!  inutrn   : Transport switch.
!
!     inutrn = 0 : neutrino transport bypassed.
!     inutrn = 1 : neutrino transport included.
!
!  ireltrns : Newtonian - GR transport toggle.
!
!     ireltrns = 0 : Newtonian neutrino transport.
!     ireltrns = 1 : general relativistic neutrino transport.
!
!  idiff    : Transport boundary conditions toggle.
!
!     idiff = 0 : dc(i,k,n) computed normally at each iteration.
!     idiff = 1 : dc(i,k,n) is the average of the current and preceding
!      timestep value.
!     idiff = 2 : dc(i,k,n) = 0 for i = inumax, dc(i,k,n) computed normally
!      at each iteration for
!      i /= inumax.
!     idiff = 3 : dc(i,k,n) = 0 for i = imax, dc(i,k,n) as in idiff = 1
!      for i /= imax
!     idiff = 4 : dc(i,k,n) = 0 for all i.
!     idiff = 5 : dc(i,k,n) is computed only on the first iteration
!     idiff = 6 : dc(i,k,n) = 0 for i = imax, computed only on the first
!      iteration for i /= imax
!
!  inumin   : inner radial zone for which neutrino transport is computed.
!  inumax   : outer radial zone for which neutrino transport is computed.
!
!  itnu(n)  : transport temperature change switch.
!
!     itnu(n) = 0 : temperature change due to n-type neutrinos bypassed.
!     itnu(n) = 1 : temperature change due to n-type neutrinos included.
!
!  jexect   : transport temperature change for zone jexect switch.
!
!     jexect = 0 : computation of the temperature is normal.
!     jexect = j : if itnu = 0, computation of temperature is normal for
!      j = jexect.
!
!  iyenu    : transport electron fraction change switch.
!
!     iyenu = 0 : electron fraction change due to n-type neutrinos bypassed.
!     iyenu = 1 : electron fraction change due to n-type neutrinos included.
!
!  jexecy   : transport electron fraction change for zone jexect switch.
!
!     jexecy = 0 : computation of the temperature is normal.
!     jexecy = j : if iyenu = 0, computation of temperature is normal for
!      j = jexecy.
!-----------------------------------------------------------------------

INTEGER                                    :: inutrn
INTEGER                                    :: ireltrns
INTEGER                                    :: idiff
INTEGER                                    :: inumin
INTEGER                                    :: inumax
INTEGER, ALLOCATABLE, DIMENSION(:)         :: itnu
INTEGER                                    :: jexect
INTEGER                                    :: iyenu
INTEGER                                    :: jexecy

!-----------------------------------------------------------------------
!
!         \\\\\ LAGRANGIAN - EULERIAN TRANSPORT SWITCHES /////
!
!-----------------------------------------------------------------------
!     lag_tr  = ye:
!     eul_tr  = ye:
!  Lagrangian transport when moving grid option on, otherwise Eulerian
!   transport. (Control)
!
!     lag_tr  = ye:
!     eul_tr  = no:
!  Lagrangian transport always. (Control)
!
!     lag_tr  = no:
!     eul_tr  = ye:
!  Eulerian transport always. (Control)
!-----------------------------------------------------------------------

CHARACTER(LEN=2)                           :: lag_tr
CHARACTER(LEN=2)                           :: eul_tr

!-----------------------------------------------------------------------
!
!        \\\\\ TEMPERATURE AND ELECTRON FRACTION UPDATES /////
!
!-----------------------------------------------------------------------

!-----------------------------------------------------------------------
!  t_chg   = ye: compute matter temperature change due to source and
!                transport
!  t_chg   = no: temperature change to source and transport bypassed
!
!  ye_chg  = ye: compute electron fraction change due to source and
!                transport
!  ye_chg  = no: electron fraction change to source and transport
!                bypassed
!-----------------------------------------------------------------------

CHARACTER(LEN=2)                           :: t_chg
CHARACTER(LEN=2)                           :: ye_chg

!-----------------------------------------------------------------------
!  Absorption and emission controls
!-----------------------------------------------------------------------
!  iaefnp   : emission and absorption of e-neutrinos and e-antineutrinos
!   on free neutrons and protons switch.
!
!     iaefnp = 0: emission and absorption of e-neutrinos and
!      e-antineutrinos on free neutrons and protons turned off.
!     iaefnp = 1: emission and absorption of e-neutrinos and
!      e-antineutrinos on free neutrons and protons turned on.
!
!  i_aeps   : emission and aborption of e-neutrinos and e-antineutrinos
!   on free neutrons and protons phase space switch.
!
!     i_aeps = 0: emission and absorption of e-neutrinos and
!      e-antineutrinos on free neutrons and protons are computed
!      assuming no recoil, or thermal motions, and approximate
!      nucleon phase space blocking.
!     iaefnp = 1: emission and absorption of e-neutrinos and
!      e-antineutrinos on free neutrons and protons are computed
!      taking into account recoil, thermal motions, and nucleon
!      phase space blocking.
!
!  rhoaefnp : density above which emission and absorption of e-neutrinos
!   and e-antineutrinos on free neutrons and protons is turned off.
!
!  iaence   : emission and absorption of e-neutrinos on nuclei as given
!   by the FFN prescription  switch
!
!     iaence = 0: emission and absorption of e-neutrinos on nuclei as
!      given by the FFN prescription turned off.
!     iaence = 1: emission and absorption of e-neutrinos on nuclei as
!      given by the FFN prescription turned on.
!
!  edmpe    : difference between the mean excited state energy of the
!   daughter nucleus and that of the parent nucleus [MeV].
!
!  iaenca   : emission and absorption of e-antineutrinos on nuclei as
!   given by the FFN prescription  switch
!
!     iaenca = 0: emission and absorption of e-antineutrinos on nuclei as
!      given by the FFN prescription turned off.
!     iaenca = 1: emission and absorption of e-antineutrinos on nuclei as
!      given by the FFN prescription turned on.
!
!  edmpa    : difference between the mean excited state energy of the
!   daughter nucleus and that of  the parent nucleus [MeV].
!
!  iaenct: emission and absorption of e-neutrinos on nuclei as given by
!   the NSE-folded tabular data for Hix et al (2003)
!
!     iaenct  = 0: emission and absorption of e-neutrinos on nuclei as
!      given by the tabular data is turned off.
!     iaenct  = 1: emission and absorption of e-neutrinos on nuclei as
!      given by the tabular data is turned on.
!
!  roaenct: density above which emission and absorption of e-neutrinos on
!   nuclei as given by Hix et al (2003) is turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: iaefnp
INTEGER                                    :: i_aeps
INTEGER                                    :: iaence
INTEGER                                    :: iaenca
INTEGER                                    :: iaenct

REAL(double)                               :: rhoaefnp
REAL(double)                               :: edmpe
REAL(double)                               :: edmpa
REAL(double)                               :: roaenct

!-----------------------------------------------------------------------
!  Weak magnetism correction for charged current interactions
!-----------------------------------------------------------------------
!  icc_wm : charged current weak magnetism switch.
!
!     icc_wm = 0: weak magnetism correction for charged current
!      interactions turned off.
!     icc_wm = 1: weak magnetism correction for charged current
!      interactions turned on.
!-----------------------------------------------------------------------

INTEGER                                    :: icc_wm

!-----------------------------------------------------------------------
!  General scattering switch
!-----------------------------------------------------------------------
!  iscat : general scattering switch.
!
!     iscat = 0 : all scattering processes turned off.
!     iscat = 1 : any scattering process is allowed to proceed if its
!      key is turned on.
!-----------------------------------------------------------------------

INTEGER                                    :: iscat

!-----------------------------------------------------------------------
!  Neutrino nucleon scattering controls
!-----------------------------------------------------------------------
!  in     : neutrino-neutron scattering switch
!
!     in = 0 : neutrino-neutron isoenergetic scattering turned off.
!     in = 1 : neutrino-neutron isoenergetic scattering turned on.
!
!  ip     : neutrino-proton scattering switch
!
!     ip = 0: neutrino-proton isoenergetic scattering turned off.
!     ip = 1: neutrino-proton isoenergetic scattering turned on.
!
!  ietann : final state blocking switch
!
!     ietann = 0 : final state blocking due to neutron or proton final
!      state occupancy in neutrino-neutron and neutrino-proton
!      isoenergetic scattering turned off.
!     ietann = 1 : final state blocking due to neutron or proton final
!      state occupancy in neutrino-neutron and neutrino-proton
!      isoenergetic scattering turned on.
!-----------------------------------------------------------------------

INTEGER                                    :: in
INTEGER                                    :: ip
INTEGER                                    :: ietann

!-----------------------------------------------------------------------
!  Neutrino nuclei scattering controls
!-----------------------------------------------------------------------
!  ihe    : neutrino-helium isoenergetic scattering switch
!
!     ihe = 0: neutrino-helium isoenergetic scattering turned off.
!     ihe = 1: neutrino-helium isoenergetic scattering turned on.
!
!  iheavy : neutrino-heavy nucleus isoenergetic scattering switch
!
!     iheavy = 0: neutrino-heavy nucleus isoenergetic scattering turned off.
!     iheavy = 1: neutrino-heavy nucleus isoenergetic scattering turned on.
!
!  iicor  : ion-ion correlation switch.
!
!     iicor = 0: neutrino-heavy nucleus isoenergetic scattering
!      correction due to ion-ion correlation effects turned off.
!     iicor = 1: neutrino-heavy nucleus isoenergetic scattering
!      correction due to ion-ion correlation effects turned on.
!-----------------------------------------------------------------------

INTEGER                                    :: ihe
INTEGER                                    :: iheavy
INTEGER                                    :: iicor

!-----------------------------------------------------------------------
!  Weak magnetism correction for neutral current interactions
!-----------------------------------------------------------------------
!  inc_wm : neutral current weak magnetism switch.
!
!     inc_wm = 0: weak magnetism correction for neutral current
!      interactions turned off.
!     inc_wm = 1: weak magnetism correction for neutral current
!      interactions turned on.
!-----------------------------------------------------------------------

INTEGER                                    :: inc_wm

!-----------------------------------------------------------------------
!  Neutrino-electron scattering controls
!-----------------------------------------------------------------------
!  nes      : neutrino-electron scattering switch.
!
!     nes = 0 : neutrino-electron scattering turned off.
!     nes = 1 : neutrino-electron scattering turned on.
!
!  rhonesmn : density below which n-neutrino-electron scattering is
!   turned off.
!
!  rhonesmx : density above which n-neutrino-electron scattering is
!   turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: nes

REAL(double)                               :: rhonesmn
REAL(double)                               :: rhonesmx

!-----------------------------------------------------------------------
!  Neutrino-positron scattering controls
!-----------------------------------------------------------------------
!  nps      : neutrino-positron scattering switch.
!
!     nes = 0 : neutrino-positron scattering turned off.
!     nes = 1 : neutrino-positron scattering turned on.
!
!  rhonpsmn : density below which n-neutrino-positron scattering is
!              turned off.
!
!  rhonpsmx : density above which n-neutrino-positron scattering is
!              turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: nps

REAL(double)                               :: rhonpsmn
REAL(double)                               :: rhonpsmx

!-----------------------------------------------------------------------
!  Neutrino-nucleon elastic scattering controls
!-----------------------------------------------------------------------
!  isctn      : neutrino-nucleon inelastic scattering switch.
!
!     isctn = 0: neutrino-nucleon inelastic scattering turned off.
!     isctn = 1: neutrino-nucleon inelastic scattering turned on.
!                Reddy et al.
!     isctn = 2: neutrino-nucleon inelastic scattering turned on.
!                Hannested & Raffelt
!
!  rhosctnemn : density below which n-neutrino-neutrino-nucleon,
!                elastic scattering is turned off.
!
!  rhosctnemx : density above which n-neutrino-neutrino-nucleon,
!                elastic scattering is turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: isctn

REAL(double)                               :: rhosctnmn
REAL(double)                               :: rhosctnmx

!-----------------------------------------------------------------------
!  Neutrino-nucleus inelastic scattering controls
!-----------------------------------------------------------------------
!  isctnA      : neutrino-nuclei inelastic scattering switch.
!
!     isctnA = 0 : neutrino-nucleus inelastic scattering turned off.
!     isctnA = 1 : neutrino-nucleus inelastic scattering turned on.
!
!  rhosctnAemn : density below which e-neutrino-neutrino-nucleus,
!   e-antineutrino-neutrino-nucleus  inelastic scattering is turned off.
!
!  rhosctnAemx : density above which e-neutrino-neutrino-nucleus,
!   e-antineutrino-neutrino-nucleus inelastic scattering is turned off.
!
!  rhosctnAtmn : density below which x-neutrino-neutrino-nucleus,
!   x-antineutrino-neutrino-nucleus  inelastic scattering is turned off.
!
!  rhosctnAtmx : density above which x-neutrino-neutrino-nucleus,
!   x-antineutrino-neutrino-nucleus scattering is turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: isctnA

REAL(double)                               :: rhosctnAemn
REAL(double)                               :: rhosctnAemx
REAL(double)                               :: rhosctnAtmn
REAL(double)                               :: rhosctnAtmx

!-----------------------------------------------------------------------
!  Electron-positron pair annihilation controls
!-----------------------------------------------------------------------
!  ipair      : electron-positron pair annihilation switch.
!
!     ipair = 0 : electron-positron pair annihilation process turned off.
!     ipair = 1 : electron-positron pair annihilation process turned on.
!
!  rhopairemn : density below which e-neutrino - e-antineutrino pair
!   production by electron-positron pair annihilation is turned off.
!
!  rhopairemx : density above which e-neutrino - e-antineutrino pair
!   production by electron-positron pair annihilation is turned off.
!
!  rhopairtmn : density below which x-neutrino - x-antineutrino pair
!   production by electron-positron pair annihilation is turned off..
!
!  rhopairtmx : density above which x-neutrino - x-antineutrino pair
!   production by electron-positron pair annihilation is turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: ipair

REAL(double)                               :: rhopairemn
REAL(double)                               :: rhopairemx
REAL(double)                               :: rhopairtmn
REAL(double)                               :: rhopairtmx

!-----------------------------------------------------------------------
!  Bremsstrahlung pair annihilation controls
!-----------------------------------------------------------------------
!  ibrem      : bremsstrahlung pair annihilation switch.
!
!     ibrem = 0 : bremsstrahlung pair annihilation process turned off.
!     ibrem = 1 : bremsstrahlung pair annihilation process turned on.
!
!  rhobrememn : density below which e-neutrino - e-antineutrino pair
!   production by bremsstrahlung pair annihilation is turned off.
!
!  rhobrememx : density above which e-neutrino - e-antineutrino pair
!   production by bremsstrahlung pair annihilation is turned off.
!-----------------------------------------------------------------------

INTEGER                                    :: ibrem

REAL(double)                               :: rhobrememn
REAL(double)                               :: rhobrememx

!-----------------------------------------------------------------------
!  NSE flashing and deflashing controls
!-----------------------------------------------------------------------
!  tnse  : temperature at which material is flashed to nse [K]
!
!  tdnse : temperature at which material is deflashed from nse [K]
!
!  deflash_pre_bounce : deflashing permitted before bounce
!-----------------------------------------------------------------------

REAL(double)                               :: tnse
REAL(double)                               :: tdnse

LOGICAL                                    :: deflash_pre_bounce

!-----------------------------------------------------------------------
!  Time window averaging of the GR pseudo-potential
!-----------------------------------------------------------------------
!  t_wndw_grpot  : time window for averaging the GR pseudo-potential
!   [micro-seconds]
!
!  t_apply_gr_ave: time after bounce at which time-window averaging
!   is to be applied [milli-seconds]
!-----------------------------------------------------------------------

REAL(double)                               :: t_wndw_grpot
REAL(double)                               :: t_apply_gr_ave

!-----------------------------------------------------------------------
!  ave_pole: pole average width
!-----------------------------------------------------------------------

INTEGER                                    :: ave_pole

END module prb_cntl_module

max_stars_repo_name: srichers/weaklib


text:
Lake Shasta is a nice destination for a wide variety of outdoor activities.  It is located about 3 hours Northern California north of Davis on Interstate 5 in  wait for it  Shasta County.  Shasta is Californias largest reservoir and its third largest body of water (after Lake Tahoe and the Salton Sea).

If youre looking for a closer lakey option, check out Lake Berryessa.  If youre in need of Shastas closer to home, theres always Shasta Drive.  Like Lake Shasta and Mount Shasta, the road was almost certainly named for wiki:wikipedia:Shasta (soft drink) Shasta soda (available at a Nugget near you!), and probably has nothing to do with the wiki:wikipedia:Shasta (tribe) Shasta Tribe that makes Northern California their home.  Right?

 See also wiki:redding:Lake Shasta on the wiki:redding Redding Wiki


text:
      subroutine difequ ( mode, xx, v )
c  from  * a practical guide to splines *  by c. de boor
calls  ppvalu(interv)
c  to be called by   c o l l o c ,  p u t i t
c  information about the differential equation is dispensed from here
c
c******  i n p u t  ******
c  mode  an integer indicating the task to be performed.
c	 = 1   initialization
c	 = 2   evaluate  de  at  xx
c	 = 3   specify the next side condition
c	 = 4   analyze the approximation
c  xx a point at which information is wanted
c
c******  o u t p u t  ******
c  v  depends on the  mode  . see comments below
c
c     parameter npiece=100,ncoef=2000
      integer mode,   i,iside,itermx,k,kpm,l,m
      real v(20),xx,   break,coef,eps,ep1,ep2,error,factor,rho,solutn
     *		      ,s2ovep,un,x,xside
      real ppvalu
c     common /approx/ break(npiece),coef(ncoef),l,kpm
      common /approx/ break(100),coef(2000),l,kpm
      common /side/ m,iside,xside(10)
      common /other/ itermx,k,rho(19)
c
c  this sample of  difequ  is for the example in chapter xv. it is a
c  nonlinear second order two point boundary value problem.
c
					go to (10,20,30,40),mode
c  initialize everything
c  i.e. set the order  m  of the dif.equ., the nondecreasing sequence
c  xside(i),i=1,...,m, of points at which side cond.s are given and
c  anything else necessary.
   10 m = 2
      xside(1) = 0.
      xside(2) = 1.
c  *** print out heading
      print 499
  499 format(37h carrier,s nonlinear perturb. problem)
      eps = .5e-2
      print 610, eps
  610 format(5h eps ,e20.10)
c  *** set constants used in formula for solution below.
      factor = (sqrt(2.) + sqrt(3.))**2
      s2ovep = sqrt(2./eps)
c  *** initial guess for newton iteration. un(x) = x*x - 1.
      l = 1
      break(1) = 0.
      do 16 i=1,kpm
   16	 coef(i) = 0.
      coef(1) = -1.
      coef(3) = 2.
      itermx = 10
					return
c
c  provide value of left side coeff.s and right side at  xx .
c  specifically, at  xx  the dif.equ. reads
c	 v(m+1)d**m + v(m)d**(m-1) + ... + v(1)d**0  =	v(m+2)
c  in terms of the quantities v(i),i=1,...,m+2, to be computed here.
   20 continue
      v(3) = eps
      v(2) = 0.
      un = ppvalu(break,coef,l,kpm,xx,0)
      v(1) = 2.*un
      v(4) = un**2 + 1.
					return
c
c  provide the	m  side conditions. these conditions are of the form
c	 v(m+1)d**m + v(m)d**(m-1) + ... + v(1)d**0  =	v(m+2)
c  in terms of the quantities v(i),i=1,...,m+2, to be specified here.
c  note that v(m+1) = 0  for customary side conditions.
   30 v(m+1) = 0.
					go to (31,32,39),iside
   31 v(2) = 1.
      v(1) = 0.
      v(4) = 0.
					go to 38
   32 v(2) = 0.
      v(1) = 1.
      v(4) = 0.
   38 iside = iside + 1
   39					return
c
c  calculate the error near the boundary layer at  1.
   40 continue
      print 640
  640 format(44h x, g(x)  and  g(x)-f(x)  at selected points)
      x = .75
      do 41 i=1,9
	 ep1 = exp(s2ovep*(1.-x))*factor
	 ep2 = exp(s2ovep*(1.+x))*factor
	 solutn = 12./(1.+ep1)**2*ep1 + 12./(1.+ep2)**2*ep2 - 1.
	 error = solutn - ppvalu(break,coef,l,kpm,x,0)
	 print 641,x,solutn,error
  641	 format(3e20.10)
   41	 x = x + .03125
					return
      end

max_stars_repo_name: ysBach/irafdocgen

max_issues_repo_name: kirxkirx/iraf

max_forks_repo_name: kirxkirx/iraf