Copying problem

General OpenMP discussion

Copying problem

Postby biernack » Fri Oct 25, 2013 4:57 am

Hi,

I am fairly new to openMP and at the moment I'm trying to parallelize a code which does some computations on the image. It's a big image, so it's being read from the file strip by strip. Then each strip is divided into squares which then are being modified by the code. After that is done, squares are copied back to strips, which are then being saved to a new file.

I parallelized the code in a way that each thread processes one square at a time. Strip is divided to 47 squares, so if I use 2 threads 24 squares are processed by one thread and another 23 by second thread. The problem only squares processed by 'master' thread are correct, the rest is not read/copied/wrote correctly; please see the attachment.

Here is the important part of the code:
Code: Select all
      do iysub=1,nysub
c get and b/g subtract strip between iy1 and iy2; pad with zeroes to width m if needed
         istat=igetsecpad(51,nx,ny,1,iy1_array(iysub),nx,m,zbufin,nx*m,bg)
         if (istat.ne.0.) then
           write(0,*) 'igetsecpad failed. code = ',istat
           stop
         endif
c         write(0,*) 'Read lines ',iy1,iy2,' of ',ny
         my=iy2_array(iysub)-iy1_array(iysub)+1
!$OMP PARALLEL     
!$OMP& PRIVATE(ix1,ix2,mx,work1,work2,zz,out,id)
!$OMP& FIRSTPRIVATE(zbufin)
!$OMP& SHARED(zbufout,nx,ny,bg,iy1_array,xmax,ymax,ntap,filt,nfit,nmask,naxlen,my,istat)
!$OMP DO SCHEDULE(DYNAMIC)
         do ixsub=1,nxsub
            if (ixsub.eq.1) then
               ix1=1-nmask
            else
               ix1=ix1+mpatch
            endif
            ix2=ix1+m-1
            if (ix2.gt.nx+nmask) then
               ix2=nx+nmask
               ix1=max(1,ix2-m+1)
            endif
            mx=ix2-ix1+1
           
c            write(0,*) 'Getting image section ',ix1,ix2,iy1,iy2
            call frombuf(nx,m,zbufin,ix1,1,mx,my,zz)
c            istat=igetsecpad(51,nx,ny,ix1,iy1,mx,my,zz,mx*my,bg)

c convolve with shmap here
c taper the array to reduce aliasing
            call taperedge(zz,mx,my,ntap)
            call shmapxpix(filt,m,nfit,zz,out,ix1,iy1_array(iysub),xmax,ymax,work1,work2)
c
c convolution done, result in out
c            write(0,*) 'Writing subsection ',ixsub,iysub,ix1,iy1
c pick out the part of out that is unaffected by aliasing (strip off border)
            call stripborder(out,zz,mx-2*nmask,my-2*nmask,m,nmask)
            call subbg(zz,(mx-2*nmask)*(my-2*nmask),-bg)

c            write(0,*) 'Write to section starting at ',ix1,iy1
!$OMP CRITICAL
            call tobuf(nx,ny,zbufout,ix1+nmask,1,mx-2*nmask,my-2*nmask,zz)
!$OMP END CRITICAL
c            call putsec(61,zz,mx-2*nmask,my-2*nmask,ix1+nmask,iy1+nmask,naxlen)
         end do
!$OMP END DO
!$OMP END PARALLEL
c         write(0,*) 'Write out lines ',iy1+nmask,iy1-nmask+my-1,' of ',nx,ny
         call putsec(61,zbufout,nx,my-2*nmask,1,iy1_array(iysub)+nmask,naxlen)
         write(0,*) 'Finished strip ',iysub,' / ',nysub,
     :       ' = pix ',iy1_array(iysub)+nmask,' - ',iy1_array(iysub)+my-nmask-1
      end do


nxsub and nysub determines number of rows and columns (47x47 here). frombuf and tobuf pull out data from the strip and put back, respectively. If there is more explanation needed, please let me know. How can I fix my problem?
I'm using gfortran on Fedora 18.

Thanks!
Pawel
Attachments
ds9.png
Part of the output picture showing the problem. (it's a star field)
ds9.png (465.35 KiB) Viewed 2178 times
biernack
 
Posts: 6
Joined: Thu Oct 24, 2013 5:53 am
Location: The Netherlands

Re: Copying problem

Postby ftinetti » Fri Oct 25, 2013 5:43 am

Hi Pawel,

Would you please send how the data are declared and the list of "intent" in, out, and inout dummies per subroutine?

Fernando.
ftinetti
 
Posts: 567
Joined: Wed Feb 10, 2010 2:44 pm

Re: Copying problem

Postby biernack » Fri Oct 25, 2013 7:56 am

Variables defined:
real kercoef(mmpoly2,msh-2)
parameter (m=512,maxnx=30000)
dimension zz(m,m),out(m,m),zbufin(m*maxnx),zbufout(m*maxnx)
dimension work1(m*m),work2(m*m)
dimension naxlen(7)
complex filt(m/2,m,mmpoly2)
data naxlen /7*0./
integer iy1_array(47), iy2_array(47),id

Subroutines:
frombuf(nx,ny,zbuf,ix1,iy1,mx,my,zz) - IN: zbuf, OUT: zz, rest are just parameters
tobuf(nx,ny,zbuf,ix1,iy1,mx,my,zz) - IN: zz, OUT: zbuf, rest are just parameters
shmapxpix(filt,m,nfit,cpix,cout,ix1,iy1,xm,ym,ftim,ctmp) - filt, m, nfit - parameters, IN: cpix, OUT: cout, rest are parameters or variables used only in this subroutine
subbg(zz,m,bg) - IN/OUT: zz, IN: m and bg
taperedge(zz,mx,my,ntap) - IN/OUT: zz, rest are just parameters
stripborder(full,stripped,mx,my,m,nmask) - IN: full, OUT: stripped, rest are just parameters

I hope that's what you needed.

Thanks,
Pawel
biernack
 
Posts: 6
Joined: Thu Oct 24, 2013 5:53 am
Location: The Netherlands

Re: Copying problem

Postby ftinetti » Fri Oct 25, 2013 9:42 am

Hi again,

Thanks for the information. Just a few comments:
1) You could use zbufin as shared, because it is not modified in the parallel region.
2) What do you mean by "are just parameters"? Are those data modified in the subroutine or not?
3) Please run your program with openmp but setting OMP_NUM_THREADS to 1

Fernando.
ftinetti
 
Posts: 567
Joined: Wed Feb 10, 2010 2:44 pm

Re: Copying problem

Postby biernack » Fri Oct 25, 2013 10:06 am

ftinetti wrote:Hi again,

Thanks for the information. Just a few comments:
1) You could use zbufin as shared, because it is not modified in the parallel region.
2) What do you mean by "are just parameters"? Are those data modified in the subroutine or not?
3) Please run your program with openmp but setting OMP_NUM_THREADS to 1

Fernando.


So as of 1) I had a version with zbufin being shared, I just wanted to remove the risk of problem with shared memory and access from to threads, as each time they read over 200k numbers.

With 2) I meant that they are either constants or integers which are not being change.

And 3) I did that before and it runs fine, gives the correct image.

Thanks for having time to look at it,
Pawel
biernack
 
Posts: 6
Joined: Thu Oct 24, 2013 5:53 am
Location: The Netherlands

Re: Copying problem

Postby ftinetti » Fri Oct 25, 2013 2:19 pm

I'm trying to reconstruct the code and data, so please check if I'm going fine:

Code: Select all
c      Data decl.         
           real kercoef(mmpoly2,msh-2)
           parameter (m=512,maxnx=30000)
           dimension zz(m,m),out(m,m),zbufin(m*maxnx),zbufout(m*maxnx)
           dimension work1(m*m),work2(m*m)
           dimension naxlen(7)
           complex filt(m/2,m,mmpoly2)
           data naxlen /7*0./
           integer iy1_array(47), iy2_array(47),id
c ...
           do iysub=1,nysub
c            get and b/g subtract strip between iy1 and iy2; pad with zeroes to width m if needed
              istat=igetsecpad(51,nx,ny,1,iy1_array(iysub),nx,m,zbufin,nx*m,bg)
              if (istat.ne.0.) then
                write(0,*) 'igetsecpad failed. code = ',istat
                stop
              endif
c              write(0,*) 'Read lines ',iy1,iy2,' of ',ny
              my=iy2_array(iysub)-iy1_array(iysub)+1

!$OMP PARALLEL
!$OMP& PRIVATE(ix1,ix2,mx,work1,work2,zz,out,id)
!$OMP& SHARED(zbufin,zbufout,nx,ny,bg,iy1_array,xmax,ymax,ntap,filt,nfit,nmask,naxlen,my,istat)
!$OMP DO SCHEDULE(DYNAMIC)
              do ixsub=1,nxsub
                 if (ixsub.eq.1) then
                    ix1=1-nmask
                 else
                    ix1=ix1+mpatch
                 endif
                 ix2=ix1+m-1
                 if (ix2.gt.nx+nmask) then
                    ix2=nx+nmask
                    ix1=max(1,ix2-m+1)
                 endif
                 mx=ix2-ix1+1

c                 write(0,*) 'Getting image section ',ix1,ix2,iy1,iy2

c                frombuf(nx,ny,zbufin,ix1,iy1,mx,my,zz)
c                            IN: zbufin, OUT: zz, rest are just parameters
                 call frombuf(nx,m,zbufin,ix1,1,mx,my,zz)

c                 istat=igetsecpad(51,nx,ny,ix1,iy1,mx,my,zz,mx*my,bg)

c convolve with shmap here
c taper the array to reduce aliasing

c                taperedge(zz,mx,my,ntap)
c                                IN/OUT: zz, rest are just parameters
                 call taperedge(zz,mx,my,ntap)

c                 shmapxpix(filt,m,nfit,cpix, out,ix1,iy1,xm,ym,ftim,ctmp)
c                                  filt, m, nfit parameters,
c                                  IN: cpix, OUT: out, rest are parameters or variables used only in this subroutine
                 call shmapxpix(filt,m,nfit,zz,out,ix1,iy1_array(iysub),xmax,ymax,work1,work2)

c
c convolution done, result in out
c                 write(0,*) 'Writing subsection ',ixsub,iysub,ix1,iy1
c pick out the part of out that is unaffected by aliasing (strip off border)

c               stripborder(full,stripped,mx,my,m,nmask) 
c                                  IN: full, OUT: stripped, rest are just parameters
                 call stripborder(out,zz,mx-2*nmask,my-2*nmask,m,nmask)

c               subbg(zz,m,bg)
c                          IN/OUT: zz, IN: m and bg
                call subbg(zz,(mx-2*nmask)*(my-2*nmask),-bg)

c                write(0,*) 'Write to section starting at ',ix1,iy1

!$OMP CRITICAL
c               tobuf(nx,ny,zbuf,ix1,iy1,mx,my,zz)
c                       IN: zz, OUT: zbufout, rest are just parameters
                 call tobuf(nx,ny,zbufout,ix1+nmask,1,mx-2*nmask,my-2*nmask,zz)
!$OMP END CRITICAL

c                 call putsec(61,zz,mx-2*nmask,my-2*nmask,ix1+nmask,iy1+nmask,naxlen)
              end do
!$OMP END DO
!$OMP END PARALLEL

c             write(0,*) 'Write out lines ',iy1+nmask,iy1-nmask+my-1,' of ',nx,ny
              call putsec(61,zbufout,nx,my-2*nmask,1,iy1_array(iysub)+nmask,naxlen)
              write(0,*) 'Finished strip ',iysub,' / ',nysub, 
         :       ' = pix ',iy1_array(iysub)+nmask,' - ',iy1_array(iysub)+my-nmask-1
           end do


I think there could be some problem in
tobuf(nx,ny,zbufout,ix1+nmask,1,mx-2*nmask,my-2*nmask,zz)

on writing to zbufout, but I can't figure it out...

Fernando.
ftinetti
 
Posts: 567
Joined: Wed Feb 10, 2010 2:44 pm

Re: Copying problem

Postby biernack » Fri Oct 25, 2013 2:30 pm

Seems OK. If you wish I could send you the code, if that would help, but I'm not entirely happy to share the whole thing publicly.

Well, as I see it is that all non-master threads (if I use more threads than 2 then only squares written by one thread are fine) have either problem with reading from zbufin or writing to zbufout. I tried to rule out the first option by making firstprivate(zbufin). This in turn means that data is not copied into zbufout. I started to wonder if there is no variable that would be crucial in tobuf and had undeclared value in non-master threads.

Thanks,
Pawel
biernack
 
Posts: 6
Joined: Thu Oct 24, 2013 5:53 am
Location: The Netherlands

Re: Copying problem

Postby ftinetti » Fri Oct 25, 2013 2:38 pm

I started to wonder if there is no variable that would be crucial in tobuf and had undeclared value in non-master threads.


I think it is a good guess, you could try including every used data (not already declared neither as private nor as shared) as shared, so that there is not any possibility that the compiler generates some private according to default rules (which, in turn, would contain garbage in non-master threads.

HTH,

Fernando.
ftinetti
 
Posts: 567
Joined: Wed Feb 10, 2010 2:44 pm

Re: Copying problem

Postby biernack » Fri Oct 25, 2013 2:51 pm

I think I might have found it. It seems like each square relies on ix1 and previous values it took (in some cases), so I need to find a way to establish it value without relying on previous square.

Thanks for the hint where to look!

Pawel
biernack
 
Posts: 6
Joined: Thu Oct 24, 2013 5:53 am
Location: The Netherlands

Re: Copying problem

Postby ftinetti » Fri Oct 25, 2013 2:54 pm

...previous values...

could be a problem... do you have some save data in the subroutine? Or, better, you may send every subroutine data declaration, just in case...

Edit: or maybe some COMMON data usage...

Fernando.
Last edited by ftinetti on Fri Oct 25, 2013 2:58 pm, edited 1 time in total.
ftinetti
 
Posts: 567
Joined: Wed Feb 10, 2010 2:44 pm

Next

Return to Using OpenMP

Who is online

Users browsing this forum: Exabot [Bot] and 6 guests