Advects a vector field.
93 int nqtot = fields[0]->GetTotPoints();
94 ASSERTL1(nConvectiveFields == inarray.num_elements(),
"Number of convective fields and Inarray are not compatible");
97 int ndim = advVel.num_elements();
98 Array<OneD, Array<OneD, NekDouble> > AdvVel (advVel.num_elements());
100 Array<OneD, Array<OneD, NekDouble> > velocity(ndim);
101 for(
int i = 0; i < ndim; ++i)
106 velocity[i] = Array<OneD, NekDouble>(nqtot,0.0);
107 fields[i]->HomogeneousBwdTrans(advVel[i],velocity[i]);
111 velocity[i] = advVel[i];
115 int nPointsTot = fields[0]->GetNpoints();
116 Array<OneD, NekDouble> grad0,grad1,grad2,wkSp;
123 nPointsTot = fields[0]->Get1DScaledTotPoints(OneDptscale);
129 for(
int i = 0; i < ndim; ++i)
131 AdvVel[i] = Array<OneD, NekDouble> (nPointsTot);
133 fields[0]->PhysInterp1DScaled(OneDptscale,velocity[i],AdvVel[i]);
138 for(
int i = 0; i < ndim; ++i)
140 AdvVel[i] = velocity[i];
144 wkSp = Array<OneD, NekDouble> (nPointsTot);
150 grad0 = Array<OneD, NekDouble> (fields[0]->GetNpoints());
151 for(
int n = 0; n < nConvectiveFields; ++n)
153 fields[0]->PhysDeriv(inarray[n],grad0);
156 Array<OneD, NekDouble> Outarray(nPointsTot);
157 fields[0]->PhysInterp1DScaled(OneDptscale,grad0,wkSp);
158 Vmath::Vmul (nPointsTot,wkSp,1,AdvVel[0],1,Outarray,1);
160 fields[0]->PhysGalerkinProjection1DScaled(OneDptscale,Outarray,outarray[n]);
164 Vmath::Vmul(nPointsTot,grad0,1,AdvVel[0],1,outarray[n],1);
169 grad0 = Array<OneD, NekDouble> (fields[0]->GetNpoints());
170 grad1 = Array<OneD, NekDouble> (fields[0]->GetNpoints());
171 for(
int n = 0; n < nConvectiveFields; ++n)
173 fields[0]->PhysDeriv(inarray[n],grad0,grad1);
177 Array<OneD, NekDouble> Outarray(nPointsTot);
178 fields[0]->PhysInterp1DScaled(OneDptscale,grad0,wkSp);
179 Vmath::Vmul (nPointsTot,wkSp,1,AdvVel[0],1,Outarray,1);
180 fields[0]->PhysInterp1DScaled(OneDptscale,grad1,wkSp);
181 Vmath::Vvtvp(nPointsTot,wkSp,1,AdvVel[1],1,Outarray,1,Outarray,1);
183 fields[0]->PhysGalerkinProjection1DScaled(OneDptscale,Outarray,outarray[n]);
187 Vmath::Vmul (nPointsTot,grad0,1,AdvVel[0],1,outarray[n],1);
188 Vmath::Vvtvp(nPointsTot,grad1,1,AdvVel[1],1,outarray[n],1,outarray[n],1);
195 Array<OneD, Array<OneD, NekDouble> > grad (ndim);
196 Array<OneD, Array<OneD, NekDouble> > gradScaled (ndim*nConvectiveFields);
197 Array<OneD, Array<OneD, NekDouble> > Outarray (nConvectiveFields);
198 for (
int i = 0; i < ndim; i++)
200 grad[i] = Array<OneD, NekDouble>(fields[0]->GetNpoints());
202 for (
int i = 0; i < ndim*nConvectiveFields; i++)
204 gradScaled[i] = Array<OneD, NekDouble>(nPointsTot);
206 for (
int i = 0; i < nConvectiveFields; i++)
208 Outarray[i] = Array<OneD, NekDouble>(nPointsTot);
211 for (
int n = 0; n < nConvectiveFields; n++)
213 fields[0]->PhysDeriv(inarray[n],grad[0],grad[1],grad[2]);
214 for (
int i = 0; i < ndim; i++)
216 fields[0]->PhysInterp1DScaled(OneDptscale,grad[i],
217 gradScaled[n*ndim+i]);
221 fields[0]->DealiasedDotProd(AdvVel,gradScaled,Outarray,
m_CoeffState);
223 for (
int n = 0; n < nConvectiveFields; n++)
225 fields[0]->PhysGalerkinProjection1DScaled(OneDptscale,
226 Outarray[n],outarray[n]);
231 Array<OneD, Array<OneD, NekDouble> > grad (ndim*nConvectiveFields);
232 Array<OneD, Array<OneD, NekDouble> > Outarray (nConvectiveFields);
233 for (
int i = 0; i < ndim*nConvectiveFields; i++)
235 grad[i] = Array<OneD, NekDouble>(nPointsTot);
237 for (
int i = 0; i < nConvectiveFields; i++)
239 Outarray[i] = Array<OneD, NekDouble>(nPointsTot);
242 for (
int n = 0; n < nConvectiveFields; n++)
244 fields[0]->PhysDeriv(inarray[n],grad[n*ndim+0],
249 fields[0]->DealiasedDotProd(AdvVel,grad,outarray,
m_CoeffState);
253 grad0 = Array<OneD, NekDouble> (fields[0]->GetNpoints());
254 grad1 = Array<OneD, NekDouble> (fields[0]->GetNpoints());
255 grad2 = Array<OneD, NekDouble> (fields[0]->GetNpoints());
256 for(
int n = 0; n < nConvectiveFields; ++n)
258 if (fields[0]->GetWaveSpace() ==
true &&
264 fields[0]->PhysDeriv(velocity[n],grad0,grad1);
268 fields[0]->HomogeneousBwdTrans(inarray[n],wkSp);
269 fields[0]->PhysDeriv(wkSp,grad0,grad1);
275 fields[0]->HomogeneousBwdTrans(outarray[n],grad2);
277 else if (fields[0]->GetWaveSpace() ==
true &&
283 fields[0]->PhysDeriv(velocity[n],grad0);
287 fields[0]->HomogeneousBwdTrans(inarray[n],wkSp);
288 fields[0]->PhysDeriv(wkSp,grad0);
293 fields[0]->HomogeneousBwdTrans(outarray[n],grad1);
297 fields[0]->HomogeneousBwdTrans(outarray[n],grad2);
301 fields[0]->PhysDeriv(inarray[n],grad0,grad1,grad2);
305 Array<OneD, NekDouble> Outarray(nPointsTot);
306 fields[0]->PhysInterp1DScaled(OneDptscale,grad0,wkSp);
307 Vmath::Vmul(nPointsTot,wkSp,1,AdvVel[0],1,Outarray,1);
309 fields[0]->PhysInterp1DScaled(OneDptscale,grad1,wkSp);
313 fields[0]->PhysInterp1DScaled(OneDptscale,grad2,wkSp);
316 fields[0]->PhysGalerkinProjection1DScaled(OneDptscale,
317 Outarray,outarray[n]);
321 Vmath::Vmul(nPointsTot,grad0,1,AdvVel[0],1,outarray[n],1);
322 Vmath::Vvtvp(nPointsTot,grad1,1,AdvVel[1],1,outarray[n],1,
324 Vmath::Vvtvp(nPointsTot,grad2,1,AdvVel[2],1,outarray[n],1,
328 if(fields[0]->GetWaveSpace() ==
true)
330 fields[0]->HomogeneousFwdTrans(outarray[n],outarray[n]);
336 ASSERTL0(
false,
"dimension unknown");
339 for(
int n = 0; n < nConvectiveFields; ++n)
bool m_homogen_dealiasing
#define ASSERTL0(condition, msg)
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
void Neg(int n, T *x, const int incx)
Negate x = -x.
MultiRegions::CoeffState m_CoeffState
MultiRegions::Direction const DirCartesianMap[]
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.