explicit discretization for diffusion term in case of 3DHomoneneous1D case
Hi All, I'd like to use explicit scheme to discretize the diffusion term in the case of 3DHomogeneous1D in continuous Galerkin projection framework. For example, I have an array u, which stores the value of first component of velocity vector in Fourier spoace. and I want to calculate the inner product (grad w, grad u), where w is the base function, grad is gradient operator. Could following code get the correct result? m_fields[0]->PhysDeriv(0,u ,u_derv[0]); // gradient in x direction m_fields[0]->PhysDeriv(1,u, u_derv[1]); // gradient in y direction for(unsigned int i=0; i<num_planes; ++i) { double beta = 2*M_PI*m_trans->GetK(i)/m_homoLen; beta *= beta; m_fields[0]->GetPlane(i)->IProductWRTDerivBase(0, &u _derv[0]+i*np, &wkCoef1[0] +i*nc); // m_fields[0]->GetPlane(i)->IProductWRTDerivBase(1, &u _derv[1]+i*np, &wkCoef1[0]+i*nc); m_fields[0]->GetPlane(i)->IProductWRTBase(tmp1 = &u[0]+i*np, &wkCoef2[0]+i*nc,MultiRegions::eGlobal); Vmath::Smul(nc, beta, &wkCoef2[0] + i*nc, 1, &wkCoef2[0]+ i*nc, 1); } Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef1,1,diff,1); Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef2,1,diff,1); // array diff stores the final diffusion value
Dear Zhicheng, After a look over your suggestion I think your approach is generally correct. However I am not sure you should be using the “MultiRegions::eGlobal” argument in IProductWRTbase. This option is stating that it will place the global coefficients into a storage format where every global degree of freedom is only stored one (hence the eGlobal). Normally we keep the coefficient space in local or elemental format. This seems to be how you are requesting the IProductWRTDerivBase to be stored. Best regards, Spencer. On 11 Jun 2015, at 03:02, Zhicheng Wang <wangzhicheng09@gmail.com<mailto:wangzhicheng09@gmail.com>> wrote: Hi All, I'd like to use explicit scheme to discretize the diffusion term in the case of 3DHomogeneous1D in continuous Galerkin projection framework. For example, I have an array u, which stores the value of first component of velocity vector in Fourier spoace. and I want to calculate the inner product (grad w, grad u), where w is the base function, grad is gradient operator. Could following code get the correct result? m_fields[0]->PhysDeriv(0,u ,u_derv[0]); // gradient in x direction m_fields[0]->PhysDeriv(1,u, u_derv[1]); // gradient in y direction for(unsigned int i=0; i<num_planes; ++i) { double beta = 2*M_PI*m_trans->GetK(i)/m_homoLen; beta *= beta; m_fields[0]->GetPlane(i)->IProductWRTDerivBase(0, &u _derv[0]+i*np, &wkCoef1[0] +i*nc); // m_fields[0]->GetPlane(i)->IProductWRTDerivBase(1, &u _derv[1]+i*np, &wkCoef1[0]+i*nc); m_fields[0]->GetPlane(i)->IProductWRTBase(tmp1 = &u[0]+i*np, &wkCoef2[0]+i*nc,MultiRegions::eGlobal); Vmath::Smul(nc, beta, &wkCoef2[0] + i*nc, 1, &wkCoef2[0]+ i*nc, 1); } Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef1,1,diff,1); Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef2,1,diff,1); // array diff stores the final diffusion value _______________________________________________ Nektar-users mailing list Nektar-users@imperial.ac.uk<mailto:Nektar-users@imperial.ac.uk> https://mailman.ic.ac.uk/mailman/listinfo/nektar-users Spencer Sherwin McLaren Racing/Royal Academy of Engineering Research Chair, Professor of Computational Fluid Mechanics, Department of Aeronautics, Imperial College London South Kensington Campus London SW7 2AZ s.sherwin@imperial.ac.uk<mailto:s.sherwin@imperial.ac.uk> +44 (0) 20 759 45052
Hi Spencer, Thank you very much for your last advice. But I have an additional problem for 3DHomogeneous1D case. My problem is how to get correct data-offset after operation m_fields[0]->*HomogeneousBwdTrans*(). As I know, the *HomogeneousBwdTrans()* operation includes *TransposeYZtoX()*,which means it is SEM ordering, so I think I should define the data-offset as follow, for(int i=0; i<num_exp; ++i) for(int j=0; j<num_planes; ++j) for(int q=0; q<nq;++q) offset = i*num_planes*nq+j*nq+q; However, when I use above procedure to get values in physical space, it gives me wrong answer. Could you please give me some advice on this problem ? zhicheng On Sun, Jun 14, 2015 at 7:10 PM Sherwin, Spencer J <s.sherwin@imperial.ac.uk> wrote:
Dear Zhicheng,
After a look over your suggestion I think your approach is generally correct. However I am not sure you should be using the “MultiRegions::eGlobal” argument in IProductWRTbase. This option is stating that it will place the global coefficients into a storage format where every global degree of freedom is only stored one (hence the eGlobal). Normally we keep the coefficient space in local or elemental format. This seems to be how you are requesting the IProductWRTDerivBase to be stored.
Best regards, Spencer.
On 11 Jun 2015, at 03:02, Zhicheng Wang <wangzhicheng09@gmail.com> wrote:
Hi All,
I'd like to use explicit scheme to discretize the diffusion term in the case of 3DHomogeneous1D in continuous Galerkin projection framework. For example, I have an array u, which stores the value of first component of velocity vector in Fourier spoace. and I want to calculate the inner product (grad w, grad u), where w is the base function, grad is gradient operator. Could following code get the correct result?
m_fields[0]->PhysDeriv(0,u ,u_derv[0]); // gradient in x direction
m_fields[0]->PhysDeriv(1,u, u_derv[1]); // gradient in y direction
for(unsigned int i=0; i<num_planes; ++i)
{
double beta = 2*M_PI*m_trans->GetK(i)/m_homoLen;
beta *= beta;
m_fields[0]->GetPlane(i)->IProductWRTDerivBase(0, &u _derv[0]+i*np, &wkCoef1[0] +i*nc); //
m_fields[0]->GetPlane(i)->IProductWRTDerivBase(1, &u _derv[1]+i*np, &wkCoef1[0]+i*nc);
m_fields[0]->GetPlane(i)->IProductWRTBase(tmp1 = &u[0]+i*np, &wkCoef2[0]+i*nc,MultiRegions::eGlobal);
Vmath::Smul(nc, beta, &wkCoef2[0] + i*nc, 1, &wkCoef2[0]+ i*nc, 1);
}
Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef1,1,diff,1);
Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef2,1,diff,1); // array diff stores the final diffusion value
_______________________________________________ Nektar-users mailing list Nektar-users@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/nektar-users
Spencer Sherwin McLaren Racing/Royal Academy of Engineering Research Chair, Professor of Computational Fluid Mechanics, Department of Aeronautics, Imperial College London South Kensington Campus London SW7 2AZ
s.sherwin@imperial.ac.uk +44 (0) 20 759 45052
Hi Zhicheng, I believe the ordering is not quite right: memory is traversed quickest by iterating over planes, then elements. However, you should probably not iterate over memory manually like this. We have routines to access the offsets of elemental storage inside the ExpList objects. It's better to use these since, for computational reasons, we might choose to reorder storage. (This is indeed what happens for mixed quad/tri meshes in 2D). If you don't care about the tensor product of quadrature inside each element, you can do something like: num_exp = m_fields[0]->GetPlane(0)->GetExpSize(); for (cnt = i = 0; i < num_planes; ++i) { for (j = 0; j < num_exp; ++j, ++cnt) { // Offset inside plane phys_offset = m_fields[0]->GetPhys_Offset(cnt); // Quadrature points in each direction nq = m_fields[0]->GetExp(cnt)->GetTotPoints(); for (k = 0; k < nq; ++k) { offset = phys_offset + k; } } } Disclaimer is that I've just written this without any kind of testing, but I believe it should work! Cheers, Dave On 15 Jun 2015, at 08:48, Zhicheng Wang <wangzhicheng09@gmail.com> wrote:
Hi Spencer,
Thank you very much for your last advice. But I have an additional problem for 3DHomogeneous1D case. My problem is how to get correct data-offset after operation m_fields[0]->HomogeneousBwdTrans(). As I know, the HomogeneousBwdTrans() operation includes TransposeYZtoX(),which means it is SEM ordering, so I think I should define the data-offset as follow,
for(int i=0; i<num_exp; ++i)
for(int j=0; j<num_planes; ++j)
for(int q=0; q<nq;++q)
offset = i*num_planes*nq+j*nq+q;
However, when I use above procedure to get values in physical space, it gives me wrong answer. Could you please give me some advice on this problem ?
zhicheng
On Sun, Jun 14, 2015 at 7:10 PM Sherwin, Spencer J <s.sherwin@imperial.ac.uk> wrote: Dear Zhicheng,
After a look over your suggestion I think your approach is generally correct. However I am not sure you should be using the “MultiRegions::eGlobal” argument in IProductWRTbase. This option is stating that it will place the global coefficients into a storage format where every global degree of freedom is only stored one (hence the eGlobal). Normally we keep the coefficient space in local or elemental format. This seems to be how you are requesting the IProductWRTDerivBase to be stored.
Best regards, Spencer.
On 11 Jun 2015, at 03:02, Zhicheng Wang <wangzhicheng09@gmail.com> wrote:
Hi All,
I'd like to use explicit scheme to discretize the diffusion term in the case of 3DHomogeneous1D in continuous Galerkin projection framework. For example, I have an array u, which stores the value of first component of velocity vector in Fourier spoace. and I want to calculate the inner product (grad w, grad u), where w is the base function, grad is gradient operator. Could following code get the correct result?
m_fields[0]->PhysDeriv(0,u ,u_derv[0]); // gradient in x direction
m_fields[0]->PhysDeriv(1,u, u_derv[1]); // gradient in y direction
for(unsigned int i=0; i<num_planes; ++i)
{
double beta = 2*M_PI*m_trans->GetK(i)/m_homoLen;
beta *= beta;
m_fields[0]->GetPlane(i)->IProductWRTDerivBase(0, &u _derv[0]+i*np, &wkCoef1[0] +i*nc); //
m_fields[0]->GetPlane(i)->IProductWRTDerivBase(1, &u _derv[1]+i*np, &wkCoef1[0]+i*nc);
m_fields[0]->GetPlane(i)->IProductWRTBase(tmp1 = &u[0]+i*np, &wkCoef2[0]+i*nc,MultiRegions::eGlobal);
Vmath::Smul(nc, beta, &wkCoef2[0] + i*nc, 1, &wkCoef2[0]+ i*nc, 1);
}
Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef1,1,diff,1);
Vmath::Vadd(NumTotCoeffs,wkCoef0,1,wkCoef2,1,diff,1); // array diff stores the final diffusion value _______________________________________________ Nektar-users mailing list Nektar-users@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/nektar-users
Spencer Sherwin McLaren Racing/Royal Academy of Engineering Research Chair, Professor of Computational Fluid Mechanics, Department of Aeronautics, Imperial College London South Kensington Campus London SW7 2AZ
s.sherwin@imperial.ac.uk +44 (0) 20 759 45052
_______________________________________________ Nektar-users mailing list Nektar-users@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/nektar-users
-- David Moxey (Research Associate) d.moxey@imperial.ac.uk | www.imperial.ac.uk/people/d.moxey Room 363, Department of Aeronautics, Imperial College London, London, SW7 2AZ, UK.
participants (3)
- 
                
                David Moxey
- 
                
                Sherwin, Spencer J
- 
                
                Zhicheng Wang