For the record, this makes much more sense written as a tensor contraction, rather than grouping $(w \cdot \nabla)$. The paper I'm referencing just writes that grouping and the equivalent summation, but skipped over writing it as a tensor contraction, which would have helped me a lot. I wrote down the equivalence here:
<https://lh3.googleusercontent.com/-BfOP_i7Vx_s/WOYKfP1BodI/AAAAAAAAUE0/0_3yExW3mGAvbRKhOElaF5rTccfwkHcwgCLcB/s1600/NonlinearNS_TensorContraction.PNG> On Thursday, April 6, 2017 at 10:41:07 AM UTC+2, Alex Zimmerman wrote: > > Correction! When I tested your idea in my code, I didn't swap the > positions of v and w, so I had w on the left and v on the right :) Now I > see that you're correct, and I'm sure that writing down the contraction by > hand would confirm it. > > So my question about the indexing of gradv is still open, but you saved me > from using that ugly loop. Thanks! > > On Thursday, April 6, 2017 at 10:33:46 AM UTC+2, Alex Zimmerman wrote: >> >> Thanks for following up. >> >> I don't think that >> double sum = v*gradz*w >> >> is equivalent to >> double sum = 0.; >> for (unsigned int i = 0; i < dim ; ++i) >> { >> for (unsigned int j = 0; j < dim; ++j) >> { >> sum += _w[j]*_gradz[i][j]*_v[i]; >> } >> } >> >> Most importantly, making the change definitely changes the behavior of my >> code; but I think this is explained by the term I wrote before which >> convects the velocity field, i.e. $[(\mathbf{w} \cdot \nabla) \mathbf{z}] >> \cdot \mathbf{v}$. The problem is how the $\nabla$ operator is grouped with >> $w$ on the left. The notation might be specific to the continuum mechanics >> community. First the dot product must be taken between the vector w and the >> gradient operator (which requires viewing the gradient operator as a >> vector), and then this result is multiplied by z, and then the dot product >> is taken between that result and v. Of course I'm doing this for the first >> time, so I could be misunderstanding something. >> >> Later on I might want to implement this differently, maybe more like some >> of the tutorials you mentioned; but I'm trying to reproduce a result from a >> paper which uses this form. It's somewhat unique because it uses Newton's >> method to linearize the nonlinear variational form. >> >> Now that I write this, I recall that the deal.II homepage News feed links >> to someone solving the stationary problem with Newton's method, step-57 >> <http://www.dealii.org/developer/doxygen/deal.II/step_57.html>. This >> didn't show up when I was glancing through the tutorials recently. I'm >> going to see if that tutorial helps :) >> >> On Wednesday, April 5, 2017 at 8:14:36 PM UTC+2, Jean-Paul Pelteret wrote: >>> >>> Hi Alex, >>>>> >>>> >>> >>>> Thanks for the clarification, Jean-Paul. >>>> >>> >>> I'm sorry to say that you shouldn't be too quick in this instance to >>> thank me... >>> >>> >>>> I'm glad to hear that my syntax is correct. >>>> >>> >>>> I am indeed using an extractor, e.g. here's a line from my code: >>>> const Tensor<2, dim> gradv = fe_values[velocity].gradient(j, quad); >>>> >>> >>> Ok, so this demonstrates something that I came to my mind only after I'd >>> posted my reply. What you're asking for on this line is the gradient of >>> the j'th shape function at the quad'th quadrature point. So not the >>> gradient of the velocity (solution) as I'd originally thought. So you'll >>> have to scrap most of what I'd said. >>> >>> I'm not too familiar with FEM implementations of fluids problems, so I'm >>> going to have to defer to other people to help you with the specifics. I >>> wouldn't want to lead you astray... But have you seen the tutorials that >>> deal with Stokes / Navier-Stokes problems, such as step-22 >>> <https://dealii.org/8.4.1/doxygen/deal.II/step_22.html#StokesProblemassemble_system> >>> and >>> 35 >>> <https://dealii.org/8.4.1/doxygen/deal.II/step_35.html#ThecodeNavierStokesProjectionassemble_advection_termcodemethodandrelated>? >>> >>> I know that they don't use this component syntax, but they might give you a >>> better idea as to how to implement your problem in general. You might also >>> find the module on handling vector valued problems >>> <https://dealii.org/8.4.1/doxygen/deal.II/group__vector__valued.html> >>> useful >>> in terms of seeing how the structure of a fluids problem might look like >>> when using / not using extractors in its implementation. >>> >>> >>>> It seems they have to do some work to map to the correct indices, i.e. >>>> with >>>> const unsigned int component_i = fe.system_to_component_index >>>> <https://www.google.com/url?q=https%3A%2F%2Fdealii.org%2Fdeveloper%2Fdoxygen%2Fdeal.II%2FclassFiniteElement.html%23a27220a135402b96c7e6eecbb04acda56&sa=D&sntz=1&usg=AFQjCNH9FPtvXF6muGtnl4nIgc1tGhyt0w> >>>> (i).first; >>>> >>> >>> This is the vector component associated with the i'th shape function. >>> I'm not sure that this is helpful for what you're wanting to implement. >>> This function and its friends are discussed in step-8 and the module on >>> handling >>> vector valued problems >>> <https://dealii.org/8.4.1/doxygen/deal.II/group__vector__valued.html> but, >>> as I said, I don't think it applies here. >>> >>> >>>> When I implement this with higher level library, FEniCS, I don't have >>>> to write the summation; but I think I have to write the summation in >>>> deal.II; hence I need to access $\partial_j z_i$, which I'm doing with the >>>> S[i][j] syntax that we discussed. I actually implement something >>>> similar to the integrand of this operator, for a given quadrature point >>>> and >>>> degree of freedom, with the lambda function: >>>> >>> auto c = []( >>>> const Tensor<1, dim> _w, >>>> const Tensor<2, dim> _gradz, >>>> const Tensor<1, dim> _v) >>>> { >>>> double sum = 0.; >>>> >>>> for (unsigned int i = 0; i < dim ; ++i) >>>> { >>>> for (unsigned int j = 0; j < dim; ++j) >>>> { >>>> >>>> >>>> sum += _w[j]*_gradz[i][j]*_v[i]; >>>> >>>> } >>>> } >>>> >>>> return sum; >>>> }; >>>> >>> >>> No, you shouldn't have to do this summation by hand (unless the tensors >>> are incompatible with one another, which is not the case here). The tensor >>> library has a pretty comprehensive set of (well-tested) contraction >>> functions implemented for it. So you'd only need to write >>> const double sum = v*_gradz*w; >>> >>> P.S. What's the best way to communicate math equations on the mailing >>>> list? I'm just typing what I think is roughly the correct LaTeX syntax, >>>> but >>>> this seems a bit clunky. >>>> >>> >>> It seems that most of us can get the gist of a LaTex-like notation, but >>> some choose to share PDFs for more elaborate equations or deviations. I >>> don't think that there's a preferred method of communicating this type of >>> information :-/ >>> >> -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.
