Advects a vector field.
95 int ndim = advVel.num_elements();
96 int nqtot = fields[0]->GetTotPoints();
97 ASSERTL1(nConvectiveFields == inarray.num_elements(),
"Number of convective fields and Inarray are not compatible");
99 Array<OneD, Array<OneD, NekDouble> > velocity(ndim);
100 for(
int i = 0; i < ndim; ++i)
104 velocity[i] = Array<OneD, NekDouble>(nqtot,0.0);
105 fields[i]->HomogeneousBwdTrans(advVel[i],velocity[i]);
109 velocity[i] = advVel[i];
113 for(
int n = 0; n < nConvectiveFields; ++n)
117 int nPointsTot = fields[0]->GetNpoints();
118 Array<OneD, NekDouble> gradV0,gradV1,gradV2, tmp, Up;
120 gradV0 = Array<OneD, NekDouble> (nPointsTot);
121 tmp = Array<OneD, NekDouble> (nPointsTot);
127 fields[0]->PhysDeriv(inarray[n],gradV0);
128 Vmath::Vmul(nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
129 Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
130 fields[0]->PhysDeriv(gradV0,tmp);
131 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
132 Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
135 gradV1 = Array<OneD, NekDouble> (nPointsTot);
136 fields[0]->PhysDeriv(inarray[n],gradV0,gradV1);
137 Vmath::Vmul (nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
138 Vmath::Vvtvp(nPointsTot,gradV1,1,velocity[1],1,outarray[n],1,outarray[n],1);
139 Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
140 Vmath::Vmul(nPointsTot,inarray[n],1,velocity[1],1,gradV1,1);
142 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
144 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
145 Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
148 gradV1 = Array<OneD, NekDouble> (nPointsTot);
149 gradV2 = Array<OneD, NekDouble> (nPointsTot);
151 fields[0]->PhysDeriv(inarray[n],gradV0,gradV1,gradV2);
157 fields[0]->DealiasedProd(velocity[0],gradV0,gradV0,
m_CoeffState);
158 fields[0]->DealiasedProd(velocity[1],gradV1,gradV1,
m_CoeffState);
159 fields[0]->DealiasedProd(velocity[2],gradV2,gradV2,
m_CoeffState);
160 Vmath::Vadd(nPointsTot,gradV0,1,gradV1,1,outarray[n],1);
161 Vmath::Vadd(nPointsTot,gradV2,1,outarray[n],1,outarray[n],1);
162 fields[0]->DealiasedProd(inarray[n],velocity[0],gradV0,
m_CoeffState);
163 fields[0]->DealiasedProd(inarray[n],velocity[1],gradV1,
m_CoeffState);
164 fields[0]->DealiasedProd(inarray[n],velocity[2],gradV2,
m_CoeffState);
166 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
168 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
170 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
171 Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
175 Up = Array<OneD, NekDouble> (nPointsTot);
178 fields[0]->HomogeneousBwdTrans(gradV0,tmp);
179 Vmath::Vmul(nPointsTot,tmp,1,velocity[0],1,outarray[n],1);
180 fields[0]->HomogeneousBwdTrans(gradV1,tmp);
181 Vmath::Vvtvp(nPointsTot,tmp,1,velocity[1],1,outarray[n],1,outarray[n],1);
182 fields[0]->HomogeneousBwdTrans(gradV2,tmp);
183 Vmath::Vvtvp(nPointsTot,tmp,1,velocity[2],1,outarray[n],1,outarray[n],1);
185 fields[0]->HomogeneousBwdTrans(inarray[n],Up);
186 Vmath::Vmul(nPointsTot,Up,1,velocity[0],1,gradV0,1);
187 Vmath::Vmul(nPointsTot,Up,1,velocity[1],1,gradV1,1);
188 Vmath::Vmul(nPointsTot,Up,1,velocity[2],1,gradV2,1);
190 fields[0]->SetWaveSpace(
false);
192 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
194 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
196 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
197 fields[0]->SetWaveSpace(
true);
200 fields[0]->HomogeneousFwdTrans(tmp,outarray[n]);
204 Vmath::Vmul(nPointsTot,gradV0,1,velocity[0],1,outarray[n],1);
205 Vmath::Vvtvp(nPointsTot,gradV1,1,velocity[1],1,outarray[n],1,outarray[n],1);
206 Vmath::Vvtvp(nPointsTot,gradV2,1,velocity[2],1,outarray[n],1,outarray[n],1);
207 Vmath::Vmul(nPointsTot,inarray[n],1,velocity[0],1,gradV0,1);
208 Vmath::Vmul(nPointsTot,inarray[n],1,velocity[1],1,gradV1,1);
209 Vmath::Vmul(nPointsTot,inarray[n],1,velocity[2],1,gradV2,1);
211 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
213 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
215 Vmath::Vadd(nPointsTot,tmp,1,outarray[n],1,outarray[n],1);
216 Vmath::Smul(nPointsTot,0.5,outarray[n],1,outarray[n],1);
220 ASSERTL0(
false,
"Dealiasing is not allowed in combination "
221 "with the Skew-Symmetric advection form for "
222 "efficiency reasons.");
226 ASSERTL0(
false,
"dimension unknown");
#define ASSERTL0(condition, msg)
MultiRegions::CoeffState m_CoeffState
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 Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
void Neg(int n, T *x, const int incx)
Negate x = -x.
MultiRegions::Direction const DirCartesianMap[]
bool m_homogen_dealiasing
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
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.