Community
cancel
Showing results for 
Search instead for 
Did you mean: 
Highlighted
Beginner
18 Views

OpenMP parameter sweep parallel

Jump to solution

I am new to OpenMP. I want to solve a stiff ODE system for a range of parameter values using parallel do loops. I use the following Pseudocode in FORTRAN given below. However, I do not know whether calling a stiff solver(as a subroutine) inside a parallel do loop is allowed or not? Also, I want to write the time series data into files with filenames such as "par_value.txt" in the subroutine before the return to the main program. Can anyone help.

    PROGRAM OPENMP_PARALLEL_STIFF
!
      USE omp_lib
!
      IMPLICIT NONE
!
      INTEGER :: I, J
      INTEGER, PARAMETER :: RTOT=10, STOT=15
      INTEGER :: TID
      INTEGER, PARAMETER :: NUM_THREADS=3
      DOUBLE PRECISION :: T_INITIAL, T_FINAL
      CALL OMP_SET_NUM_THREADS(NUM_THREADS)
      CALL CPU_TIME(T_INITIAL)
!
     PRINT*, "TIME INITIAL ",T_INITIAL
!
!$OMP PARALLEL DO PRIVATE(I,J)
      DO I=1,RTOT
        DO J=1,STOT
        TID=OMP_GET_THREAD_NUM()
         CALL STIFF_DRIVER(TID,I,J,RTOT,STOT)
        END DO
      END DO
!$OMP END PARALLEL DO
!
      CALL CPU_TIME(T_FINAL)
!
     PRINT*, "TIME FINAL ",T_FINAL
!
     PRINT*, "TIME ELAPSED ",(T_FINAL-T_INITIAL)/NUM_THREADS
    END PROGRAM OPENMP_PARALLEL_STIFF

    SUBROUTINE STIFF_DRIVER(TID,II,JJ,RTOT,STOT)
!
      USE USEFUL_PARAMETERS_N_FUNC
      USE DVODE_F90_M

!     -----------------------------------------------------------------
!     Type declarations:
!     -----------------------------------------------------------------
      IMPLICIT NONE
!     -----------------------------------------------------------------
!     Number of odes and number of event functions for this problem:
!     -----------------------------------------------------------------
      INTEGER :: SERIAL_NUMBER, TID
      INTEGER :: II, JJ, RTOT, STOT, IND
      INTEGER :: J, NTOUT
      INTEGER :: ITASK, ISTATE, ISTATS, I
!     -----------------------------------------------------------------
!     parameters : declaration
!     -----------------------------------------------------------------
      DOUBLE PRECISION, PARAMETER :: s0=0.450D0, dr=1.0D-4, ds=1.0D-2
 !     -----------------------------------------------------------------
      DOUBLE PRECISION, DIMENSION(NEQ) :: Y, YOUT
      DOUBLE PRECISION :: ATOL, RTOL, RSTATS, T, TOUT, EPS, TFINAL, DELTAT
      DIMENSION :: RSTATS(22), ISTATS(31)
      DOUBLE PRECISION :: bb, cc, ba, ba1, eta
      CHARACTER(len=45) :: filename
!       ---------------------------------------------------------------
!
      TYPE (VODE_OPTS) :: OPTIONS
!
!     -----------------------------------------------------------------
      SERIAL_NUMBER=11+II+(JJ-1)*RTOT
!     IND=TID+11+II+(JJ-1)*RTOT
      WRITE (*,12)SERIAL_NUMBER,TID
12    FORMAT ("SL. NO. ",I3," THREAD NO.",I3)
!     -----------------------------------------------------------------
      r=(II-1)*dr
      s=s0+JJ*ds
!     -----------------------------------------------------------------

     EPS = 1.0D-9

!    ------------------------------------------------------------------
!     Open the output file and the plot files:
!     -----------------------------------------------------------------
      WRITE (filename,93)r,s
93    FORMAT ("r_",f6.4,"_s_",f4.2,".txt")
      OPEN (UNIT=NEWUNIT(IND),FILE=filename,STATUS='UNKNOWN')
!     -----------------------------------------------------------------
      EPS = 1D-9
!     -----------------------------------------------------------------
!     Parameters for the stiff ODE system
!     -----------------------------------------------------------------
       noise = 0.0010D0
       pi = 4.0D0*atan(1.0D0)
       q0 = 0.60D0;    v = 3.0D0
       Va = 20.0D-4;  Vs = 1.0D-1
       e1 = 1.0D-1;   e2 = 1.10D-5; e3 = 2.3D-3; e4=3.0D-4
      del = 1.7D-4;   mu = 5.9D-4
       al = 1.70D-4;  be = 8.9D-4;  ga = 2.5D-1
!        ---------------------------------------------------------------
!      Surrounding
!      ---------------------------------------------------------------
          e1s = e1/s;    e2s = e2/(s**2);   e3s = e3/s;    e4s = e4/s
         dels = del*s;    rs = r*s
!      --------------------------------------------------------------
!      discrete subsystems
!      --------------------------------------------------------------
          e1v = e1/v;     e2v = e2/(v**2);   e3v = e3/v;    e4v = e4/v
         delv = del*v;     rv = r*v
!      ---------------------------------------------------------------
!       Assign values for q
!       ---------------------------------------------------------------
        CALL RANDOM_SEED
         DO I = 1, n
          CALL RANDOM_NUMBER(HARVEST=bb)
          IF(bb <= 0.0D0)THEN
             bb = 0.00010D0
          END IF
          CALL RANDOM_NUMBER(HARVEST=cc)
           ba = -2.0D0*DLOG(bb)
          ba1 = DCOS(2.0D0*pi*cc)
          eta = noise*DSQRT(ba)*ba1
!       ---------------------------------------------------------------
         q(I) = q0 + eta
!       ---------------------------------------------------------------
         END DO
!     -----------------------------------------------------------------
!     Set the integration parameters:
!     -----------------------------------------------------------------
      T = 0.0D0
      TFINAL = 200.0D0
      DELTAT = 0.10D0
      NTOUT = INT(TFINAL/DELTAT)
      RTOL = EPS
      ATOL = EPS
      ITASK = 1
      ISTATE = 1
!     -----------------------------------------------------------------
!     Set the initial conditions:
!     -----------------------------------------------------------------
      CALL Y_INITIAL(NEQ,Y)
!     -----------------------------------------------------------------
!     Set the VODE_F90 options:
!     -----------------------------------------------------------------
      OPTIONS = SET_OPTS(DENSE_J=.TRUE.,USER_SUPPLIED_JACOBIAN=.FALSE., &
        RELERR=RTOL,ABSERR=ATOL,MXSTEP=100000)
!     -----------------------------------------------------------------
!     Perform the integration:
!     -----------------------------------------------------------------
      DO I=1,NTOUT
!     -----------------------------------------------------------------
      TOUT = (I-1)*DELTAT
!     -----------------------------------------------------------------
      CALL DVODE_F90(F_FUNC,NEQ,Y,T,TOUT,ITASK,ISTATE,OPTIONS)
!     -----------------------------------------------------------------
!     Stop the integration if an error occurred:
!     -----------------------------------------------------------------
      IF (ISTATE<0) THEN
        WRITE (*,*)"ISTATE ", ISTATE
        STOP
      END IF
!     -----------------------------------------------------------------
!     WRITE DATA TO FILE
!     -----------------------------------------------------------------
        WRITE (IND,*) TOUT,T, x_mean, Y(NEQ-2)
       END DO
!     -----------------------------------------------------------------
      CLOSE(UNIT=IND)


    END SUBROUTINE STIFF_DRIVER

 

I also used the following function in the modules USEFUL_PARAMETERS_N_FUNC to obtain unit number for opening files in the called subroutine to write the time series data.

      INTEGER FUNCTION NEWUNIT(UNIT)
      INTEGER, INTENT(OUT), OPTIONAL :: UNIT
!     local
      INTEGER, PARAMETER :: LUN_MIN=3011, LUN_MAX=4000
      LOGICAL :: open_ed
      INTEGER :: l
!     begin
      NEWUNIT=-1
        DO l=LUN_MIN,LUN_MAX
          INQUIRE(UNIT=l,OPENED=open_ed)
            IF(.NOT.open_ed) THEN
               NEWUNIT=l
               EXIT
            END IF
        END DO
        IF (PRESENT(UNIT)) UNIT=NEWUNIT
      END FUNCTION NEWUNIT

However, I get on compiling and running the code

Fortran runtime error: File already opened in another unit:

 A storage allocation error occurred.                                            
 The error occurred at location I1.                                              
In the above message, I1 =        280
LEVEL = 2 in XERRDV. Stopping.

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

Backtrace for this error:

0 Kudos

Accepted Solutions
Highlighted
18 Views

TID            0  filename r_0.0002_s_1.90.txt                                                        
TID            1  filename r_0.0002_s_1.90.txt                                                        
 

Two different threads are using the same file name. Did you intend to have the TID as part of the file name?

Jim Dempsey

 

View solution in original post

0 Kudos
6 Replies
Highlighted
18 Views

!$OMP PARALLEL DO PRIVATE(I,J,TID)

Prior to OPEN, print the TID and filename, following OPEN print omp_get_thread_num() and IND . Assure values are good.

Jim Dempsey

0 Kudos
Highlighted
Beginner
18 Views

Dear Jim,

I tried your suggestion and got the following on printing the OMP_GET_THREAD_NUM() and IND were also giving acceptable numbers untill the following errors occurred.

.

.

.

 TID            2  filename r_0.0069_s_1.63.txt                                                        
 TID            2  filename r_0.0069_s_1.64.txt                                                        
 TID            2  filename r_0.0069_s_1.65.txt                                                        
 TID            1  filename r_0.0036_s_1.66.txt                                                        
 TID            0  filename r_0.0002_s_1.89.txt                                                        
 TID            0  filename r_0.0002_s_1.90.txt                                                        
 TID            1  filename r_0.0002_s_1.90.txt                                                        
At line 81 of file openmp_parallel_stiff.f90 (unit = 3012)
Fortran runtime error: File already opened in another unit
 TID            0  filename r_0.0002_s_1.91.txt                                                        

Program aborted. Backtrace:
[Thread 0x7ffff6940700 (LWP 31565) exited]
[Thread 0x7ffff7fd1780 (LWP 31561) exited]
[Inferior 1 (process 31561) exited with code 02]
(gdb) step
The program is not being run.
(gdb) q

I fail to understand why in spite of using unit = base_value + TID. I get Fortran runtime error: File already opened in another unit.

When I declare TID as private should I not get distinct unit numbers using the expression above.

I  also tried RECURSIVE SUBROUTINE STIFF_DRIVER(TID,II,JJ,RTOT,STOT) but that is not helping either.

Best,

Ed

0 Kudos
Highlighted
19 Views

TID            0  filename r_0.0002_s_1.90.txt                                                        
TID            1  filename r_0.0002_s_1.90.txt                                                        
 

Two different threads are using the same file name. Did you intend to have the TID as part of the file name?

Jim Dempsey

 

View solution in original post

0 Kudos
Highlighted
Beginner
18 Views

Dear Jim,

Thank you for pointing it out the problem with the filename. So is it the filename that is creating the problem and not the UNIT number?

No I do not intend to use the TID as part of the filename. I was just printing the TID along with the filename on the screen just to monitor the working thread. It is not a part of the filename. Sorry for the confusion that the priniing format created.  However, I do not understand why filenames are same for two threads? In other words, why combination of  r and s are not unique being dependent on II and JJ (even though I,J are private in the omp parallel do directive)?

Any suggestions to get around the problem.

Best,

Ed

 

0 Kudos
Highlighted
18 Views

Nothing stands out. The only thing to suggest is to

1) run the code with 1 thread to verify it is or is not a threading issue
2) add additional diagnostic PRINT *,... (consider wrapping these in a critical section such that printouts do not intermingle)

Jim Dempsey

0 Kudos
Highlighted
Beginner
18 Views

Dear Jim,

The problem was with the format statement for r and s as you hinted when filenames were same for two threads. Sorry for the previous post. Shoudlhave been more carefull about the FORMAT specifiers.

The code is running fine now.

Thanks a lot.

Best,

Ed

 

0 Kudos