[Nek5000-users] opdiv and multd

[nek5000-users at lists.mcs.anl.gov]

dssum() is just a wrapper for call gs_op() which takes only the input field from the dssum call. But I agree in general it's a good idea to call all subroutines with the correct arguments.

-Stefan





----- Original Message -----
From: nek5000-users at lists.mcs.anl.gov
To: nek5000-users at lists.mcs.anl.gov
Sent: Tuesday, December 14, 2010 6:50:15 PM
Subject: Re: [Nek5000-users] opdiv and multd

Sorry i have to do a correction.

A user point it out that the subroutine is defined as dssum(u,nx,ny,nz) 
in dssum.f .
Thus it should be called  with 4 arguments :

      call dssum(RES(1),nx1,ny1,nz1)

This did not give me problems and it seems to me that nx,ny,nz are not 
used.
But it is better to be aware of it since i do not know how the compiler 
will read the wrong call.

Sorry about that.
francesco




On 12/08/2010 04:02 PM, nek5000-users at lists.mcs.anl.gov wrote:
>
> Thanks Francesco!
>
> Paul
>
>
> On Wed, 8 Dec 2010, nek5000-users at lists.mcs.anl.gov wrote:
>
>> Thank you Paul,
>> that solved the problem.
>>
>> The correct result is obtained by :
>>      call multd (RES(1),vx,rym2,sym2,tym2,2,iflg)
>>      call dssum(RES(1))
>>      call col2(RES(1),binvm1,lt)
>>
>> francesco
>>
>>
>>
>> On 12/08/2010 02:17 PM, nek5000-users at lists.mcs.anl.gov wrote:
>>>
>>> Hi Francesco,
>>>
>>> Several of the routines internal to nek are the so-called weak
>>> derivatives, implying that they somehow involve the mass matrix
>>> or the transpose of the operator, so it takes a bit of investigation
>>> to understand precisely what each operator is doing.  (Note that the
>>> functionality also differs slightly between PN-PN and PN-PN-2 
>>> formulations.)
>>>
>>> As an example, axhelm does not multiply by the (negative) Laplacian.
>>> Instead, it creates (for h1==1, h2==0):
>>>
>>>    w = A * u
>>>
>>>      = grad^T B grad u
>>>
>>> where B is the mass matrix.
>>>
>>> Similarly, opdiv is most likely returning   d = B*grad u (or perhaps
>>> the negative of this).     Moreover, d is returned on the pressure
>>> mesh (mesh 2), which is the same as the velocity mesh (mesh 1) for 
>>> the PN-PN formulation that you use.
>>>
>>> In Nek, the mass matrix is diagonal and the mesh 1 inverse is stored
>>> in binvm1.  It seems probable that you'll get the correct result
>>> after
>>>
>>>       call col2(d,binvm1,n)
>>>
>>> where d holds the output of opdiv and n is the number of points on
>>> mesh 1, assuming you're using PN-Pn.
>>>
>>> Paul
>>>
>>>
>>>
>>> On Wed, 8 Dec 2010, nek5000-users at lists.mcs.anl.gov wrote:
>>>
>>>> Hi all,
>>>>
>>>> I am having some problems using the subroutines opdiv and multd.
>>>> For some reasons they do not give me correct results.
>>>> Do they need any special treatment or precautions?
>>>>
>>>> I was trying to compute the divergence of a 3D field, and i 
>>>> realized that in my case the multd subroutine used
>>>> by opdiv does not gives good results, while the same derivative 
>>>> computed from gram1 gives no problems.
>>>>
>>>> For example, I overwrite an analytical field :
>>>>
>>>>         ...
>>>>         parameter(lt=lx1*ly1*lz1*lelt)
>>>>         real       Res(lt)
>>>>         ...
>>>>                    vx(i,j,k,e)=   1.0*sin(x)*cos(y)*cos(z)
>>>>         ...
>>>>             call print_line(vx)
>>>>         ...
>>>>         iflg = 1
>>>>         call multd (RES(1),vx,rym2,sym2,tym2,2,iflg)
>>>>            call print_line(Res(1) )
>>>>
>>>>        call gradm1(work1,Res(1),work2,vx)
>>>>            call print_line(Res(1) )
>>>>         ...
>>>>         if (lx2.ne.lx1) then
>>>>         ...
>>>>         STOP
>>>>         endif
>>>>
>>>>
>>>> and  print output a (xr,:,zr) the result  (performed by print_line) .
>>>> In vxi.eps i compare the result (lines) with the analytical 
>>>> solution( symbols). The result of multd is not correct while all
>>>> the others match perfectly the analytical solution.
>>>> In vxy.eps i plot only the result of multd. Note that on y i have 
>>>> 12 elements the same number of peaks of the solution.
>>>>
>>>> Does anybody have an idea on how to fix this?
>>>>
>>>> Thank you very much for any help or suggestion
>>>>
>>>> francesco
>>>>
>>>>
>>>>
>>>>
>>> _______________________________________________
>>> Nek5000-users mailing list
>>> Nek5000-users at lists.mcs.anl.gov
>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>
>> _______________________________________________
>> Nek5000-users mailing list
>> Nek5000-users at lists.mcs.anl.gov
>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>

_______________________________________________
Nek5000-users mailing list
Nek5000-users at lists.mcs.anl.gov
https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users