[Q-e-developers] Loops on real-space grid

[Oliviero Andreussi]

Thanks, indeed I overlooked that check in add_efield.f90. I will patch compute_dip.f90. For add_monopole.f90, the xy plane is not looked at explicitly, as the correction is uniform in the xy plane. Still, this means that the correction is also added to elements of the vector which are unphysical if nrx1!=nr1, which may be a bug (but it is likely not to be one, as I guess that potential is used to multiply the charge density, which is zeroed in those points). I would easily patch this file as well, but if the developers of the file have any objection, please let me know in the next couple of days. 

Concerning the more appropriate way to design the loops, I also like the second one better, shall I modify all the instances in Modules (which are in compute_dipole.f90, qmmm.f90, generate_functions.f90, and fd_gradient.f90)? The two latter subroutines are ‘mine' and I do not have any objection to change them.

Best,

Oliviero


> On 29 May 2017, at 16:36, Paolo Giannozzi <p.giannozzi at gmail.com> wrote:
> 
> On Mon, May 29, 2017 at 2:34 PM, Oliviero Andreussi
> <oliviero.andreussi at usi.ch> wrote:
> 
>> Having said this, I double checked inside PW/src and the same patch could be done
>> for PW/src/compute_dip.f90 and for PW/src/add_efield.f90
> 
> add_efield.f90 seems fine to me, thanks to this check:
> 
>     ! ... do not include points outside the physical range
> 
>     IF ( i >= dfftp%nr1 .OR. j >= dfftp%nr2 .OR. k >= dfftp%nr3 ) CYCLE
> 
> compute_dip.f90 accesses the charge density in nonphysical points and
> sums over them, but, as long as the charge density is properly set to
> 0 in nonphysical points, this should not lead to wrong results. I
> strongly prefer that nonphysical points are never accessed, though.
> 
> By the way: this problem of FFT dimensions has surfaced more than once
> in the past; it is almost forgotten since the demise of AIX machines
> with ESSL, because the trick of increasing first dimension by one is
> not implemented for FFTW and Intel FFT. Not sure whether such trick
> may still be useful in the future or whether it is a legacy of a
> distant past.
> 
>> while add_monofield already had the test in place [...]
> 
>>     ! ... do not include points outside the physical range
>> 
>>     IF ( k >= dfftp%nr3 ) CYCLE
> 
> I don't think this is sufficient, though: one should check the index
> in xy plane as well.
> 
>> I could easily patch PW/src/compute_dip.f90 (which is, btw, very very similar
>> to Modules/compute_dipole.f90, is there really a need for a duplicate?)
> 
> I have tried, and given up, more than once to merge these routines:
> they look the same but there are subtle differences.
> 
>> In particular, it seems to me that while the loops in Modules have this general preliminary setup
>> [...]
>> the same loops in PW/src have this slightly different one
>> [...]
>> I understand that the two should be equivalent, provided that all the components of the fft
>> descriptor are initialised correctly with and without MPI, so shall we only use one of the two
>> approaches? Which?
> 
> I have a slight preference for the latter because it is the same for
> MPI and non-MPI case. I would like to see #ifdef __MPI confined to
> what is really different in serial and in parallel execution
> 
> Paolo
> 
> -- 
> Paolo Giannozzi, Dip. Scienze Matematiche Informatiche e Fisiche,
> Univ. Udine, via delle Scienze 208, 33100 Udine, Italy
> Phone +39-0432-558216, fax +39-0432-558222
> _______________________________________________
> Q-e-developers mailing list
> Q-e-developers at qe-forge.org
> http://qe-forge.org/mailman/listinfo/q-e-developers