- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
hi below is the code I rewrote(with do loops) and further down(with where statements) is the original code. The code seems to be producing differnt values for WORK arrays over time. How to resolve this issue??
do j=1,ny_block
do i=1,nx_block
if(LMASK(i,j) .eq. .true. )then
WORK1(i,j,kk) = KAPPA_THIC(i,j,kbt,k,bid) &
* SLX(i,j,kk,kbt,k,bid) * dz(k)
WORK2(i,j,kk) = c2 * dzwr(k) * ( WORK1(i,j,kk) &
- KAPPA_THIC(i,j,ktp,k+1,bid) * SLX(i,j,kk,ktp,k+1,bid) &
* dz(k+1) )
WORK2_NEXT(i,j) = c2 * ( &
KAPPA_THIC(i,j,ktp,k+1,bid) * SLX(i,j,kk,ktp,k+1,bid) - &
KAPPA_THIC(i,j,kbt,k+1,bid) * SLX(i,j,kk,kbt,k+1,bid) )
WORK3(i,j,kk) = KAPPA_THIC(i,j,kbt,k,bid) &
* SLY(i,j,kk,kbt,k,bid) * dz(k)
WORK4(i,j,kk) = c2 * dzwr(k) * ( WORK3(i,j,kk) &
- KAPPA_THIC(i,j,ktp,k+1,bid) * SLY(i,j,kk,ktp,k+1,bid) &
* dz(k+1) )
WORK4_NEXT(i,j) = c2 * ( &
KAPPA_THIC(i,j,ktp,k+1,bid) * SLY(i,j,kk,ktp,k+1,bid) - &
KAPPA_THIC(i,j,kbt,k+1,bid) * SLY(i,j,kk,kbt,k+1,bid) )
endif
enddo
enddo
where ( LMASK )
WORK1(:,:,kk) = KAPPA_THIC(:,:,kbt,k,bid) &
* SLX(:,:,kk,kbt,k,bid) * dz(k)
WORK2(:,:,kk) = c2 * dzwr(k) * ( WORK1(:,:,kk) &
- KAPPA_THIC(:,:,ktp,k+1,bid) * SLX(:,:,kk,ktp,k+1,bid) &
* dz(k+1) )
WORK2_NEXT = c2 * ( &
KAPPA_THIC(:,:,ktp,k+1,bid) * SLX(:,:,kk,ktp,k+1,bid) - &
KAPPA_THIC(:,:,kbt,k+1,bid) * SLX(:,:,kk,kbt,k+1,bid) )
WORK3(:,:,kk) = KAPPA_THIC(:,:,kbt,k,bid) &
* SLY(:,:,kk,kbt,k,bid) * dz(k)
WORK4(:,:,kk) = c2 * dzwr(k) * ( WORK3(:,:,kk) &
- KAPPA_THIC(:,:,ktp,k+1,bid) * SLY(:,:,kk,ktp,k+1,bid) &
* dz(k+1) )
WORK4_NEXT = c2 * ( &
KAPPA_THIC(:,:,ktp,k+1,bid) * SLY(:,:,kk,ktp,k+1,bid) - &
KAPPA_THIC(:,:,kbt,k+1,bid) * SLY(:,:,kk,kbt,k+1,bid) )
endwhere
Are they equivalent??
Coz I am having issues of correctness here.
Link Copied
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
I would say they aren't equivalent, as the explicit DO loops forbid the compiler from switching loop nest unless it can "prove" correctness, while the multiple-rank syntax gives the compiler less guidance, and gets into insufficiently tested territory.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
I am so sorry the second one is the original code and first is changed code.
at the end of loops assuming inputs are the same do they produce same value for work arrays???
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
At the end I find work arrays are producing different values. I dont know how to ensure correctness while rewriting my code.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Try this:
where ( LMASK )
WORK1(:,:,kk) = KAPPA_THIC(:,:,kbt,k,bid) &
* SLX(:,:,kk,kbt,k,bid) * dz(k)
WORK3(:,:,kk) = KAPPA_THIC(:,:,kbt,k,bid) &
* SLY(:,:,kk,kbt,k,bid) * dz(k)
end where
where ( LMASK )
WORK2(:,:,kk) = c2 * dzwr(k) * ( WORK1(:,:,kk) &
- KAPPA_THIC(:,:,ktp,k+1,bid) * SLX(:,:,kk,ktp,k+1,bid) &
* dz(k+1) )
WORK2_NEXT = c2 * ( &
KAPPA_THIC(:,:,ktp,k+1,bid) * SLX(:,:,kk,ktp,k+1,bid) - &
KAPPA_THIC(:,:,kbt,k+1,bid) * SLX(:,:,kk,kbt,k+1,bid) )
WORK4(:,:,kk) = c2 * dzwr(k) * ( WORK3(:,:,kk) &
- KAPPA_THIC(:,:,ktp,k+1,bid) * SLY(:,:,kk,ktp,k+1,bid) &
* dz(k+1) )
WORK4_NEXT = c2 * ( &
KAPPA_THIC(:,:,ktp,k+1,bid) * SLY(:,:,kk,ktp,k+1,bid) - &
KAPPA_THIC(:,:,kbt,k+1,bid) * SLY(:,:,kk,kbt,k+1,bid) )
end where
aketh t. wrote:
At the end I find work arrays are producing different values. I dont know how to ensure correctness while rewriting my code.
You can create a minimal working example using your multi-dimensional arrays and try it out both ways, with DO loop and WHERE construct.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
err I am sorry.
I will change the main thread contents.
The first code is the changed second is the original. I want do loops in my code for openmp support.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Apart from the other issues, you may have created a subtle bug with the seemingly benign test in the line
if(LMASK(i,j) .eq. .true. )then
Presumably, LMASK is of type LOGICAL.
The issue is covered by Doctor Fortran in "To .EQV. or to .NEQV., that is the question", or "It's only LOGICAL", see https://software.intel.com/en-us/forums/topic/275071#comment-1548435 .
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
wont
if(LMASK(i,j))then
.
.
.
endif
be sufficient
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Replying to #8:
Yes; sufficient, but also LOGICAL!
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
the 2 versions may have the expressions performed in different orders, and hence you can get differing results. DO loops are optimization/vectorization candidates.
you may get closer to identical results if you use -O0
If you insist on using optimization, you can look at the optimization report via -qopt-report5 and you'll see differences in how these blocks are treated.
I can quite easily see where the results would vary slightly at optimization.
Here is a good read on this topic: https://software.intel.com/sites/default/files/managed/f1/35/fp-control-2013-03.pdf
ron
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
aketh,
The point you are missing from mecej4's post is: We know not what your LMASK is, and equally important, how you set the values into LMASK.
In Fortran (default settings) .true. is indicated by the least significant bit in the variable being set to 1. .false. is indicated by the least significant bit of the variable being 0.
Therefore, for logicals, use .EQV. and .NEQV., do not use == as when the value held in the logical variable is numerically 1, the variable may be logically .true. but not equal to logical .true..
IF(someLogical == .true. ) xxx
May (may) require the integer equivalent value within your someLogical to be -1 in order to evaluate to .true.
Whereas
IF(someLogical .EQV. .true.) xxx
Makes the decision based only on the lsb of what is in your someLogical.
As long as you do not set a logical variable with an integer value, then == and .EQV. are (likely) interchangeable..
The link that mecj4 gave you is worth reading closely.
If by chance (or intention) LMASK contains /=.true. values that are .true. then use:
if(LMASK(i,j))then
As this will evaluate only the lsb of what is inside LMASK(i,j)
Jim Dempsey
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi,
To me variation in data is okay, but what I would like is if the data are EQUIVALENT.
1.222222 and 1.222223 are okay.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
LMASK is of type logical
logical :: LMASK(x,y)
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
aketh t. wrote:
Hi,
To me variation in data is okay, but what I would like is if the data are EQUIVALENT.
1.222222 and 1.222223 are okay.
Assuming the left number is in A and the right number is in B
The expression A*EPSILON(A) will yield a number that is approximately equivalent to the value of the least significant bit in the mantissa of A. You may need (want) to use ABS(A)*EPSILON(A)
DELTA = ABS(A)*EPSILON(A)
IF(((A-DELTA).LE.B).AND.((A+DELTA).GE.B)) ...
Caution, you may need a DELTA that is larger that obtained above. It is your responsibility to determine a suitable value.
You can also use the intrinsic function NEAREST
IF((NEAREST(A,-1.0) .LE. B).AND.(NEAREST(A,1.0).GE.B)) ...
Again, convergence routines may require a larger degree of uncertainty (or certainty).
Jim Dempsey
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
hi all,
I managed to solve my issue with all your help. thanks you all.
Very often it happens we solve our problems and simply dont have the courtesy to thank you all once done.
therefore here is my special thanks to you all, without whom I would never be able go along with my efforts at parallelization and vectorization of CESM climate code.
I realize how you cater to the a huge community, which must be a thankless job.
- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page