dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

I'm not sure this will work. If you attach 6 vertices to a triangle by

  <triangle index="0" affine="false" v0="0" v1="1" v2="2" v3="4" v4="5" v5="6"/>

then all sorts of things will break (I imagine). A triangle always has
three vertices.

The geometry of the triangles is separate.

Maybe we could just add extra data which could be "control points" for
the cell facets? For P2 it would be the edge/face midpoints.

-- 
Anders


On Tue, Aug 19, 2008 at 02:21:41PM -0400, Shawn Walker wrote:
> Ok, here is what I suggest.  If anyone sees anything wrong with this  
> format, please say something!
>
> =============== THE OLD WAY =============================
>
> <?xml version="1.0" encoding="UTF-8"?>
>
> <dolfin xmlns:dolfin="http://www.fenics.org/dolfin/";>
>   <mesh celltype="triangle" dim="2">
>     <vertices size="4">
>       <vertex index="0" x="0.0" y="0.0"/>
>       <vertex index="1" x="1.0" y="0.0"/>
>       <vertex index="2" x="0.0" y="1.0"/>
>       <vertex index="3" x="-1.0" y="0.0"/>
>     </vertices>
>     <cells size="2">
>       <triangle index="0" v0="0" v1="1" v2="2"/>
>       <triangle index="1" v0="3" v1="0" v2="2"/>
>     </cells>
>   </mesh>
> </dolfin>
>
> =============== THE NEW WAY =============================
>
> # for an affinely mapped mesh; note the new variable `map_type'.
>
> <?xml version="1.0" encoding="UTF-8"?>
>
> <dolfin xmlns:dolfin="http://www.fenics.org/dolfin/";>
>   <mesh celltype="triangle" dim="2" map_type="P1">
>     <vertices size="4">
>       <vertex index="0" x="0.0" y="0.0"/>
>       <vertex index="1" x="1.0" y="0.0"/>
>       <vertex index="2" x="0.0" y="1.0"/>
>       <vertex index="3" x="-1.0" y="0.0"/>
>     </vertices>
>     <cells size="2">
>       <triangle index="0" v0="0" v1="1" v2="2"/>
>       <triangle index="1" v0="3" v1="0" v2="2"/>
>     </cells>
>   </mesh>
> </dolfin>
>
> # for a piecewise quadratic mesh; here map_type="P2".
>
> <?xml version="1.0" encoding="UTF-8"?>
>
> <dolfin xmlns:dolfin="http://www.fenics.org/dolfin/";>
>   <mesh celltype="triangle" dim="2" map_type="P2">
>     <vertices size="9">
>       <vertex index="0" x="0.0"  y="0.0"/>
>       <vertex index="1" x="1.0"  y="0.0"/>
>       <vertex index="2" x="0.0"  y="1.0"/>
>       <vertex index="3" x="-1.0" y="0.0"/>
>
>       <vertex index="4" x="0.5"  y="0.05"/>
>       <vertex index="5" x="0.6"  y="0.6"/>
>       <vertex index="6" x="0.0"  y="0.5"/>
>       <vertex index="7" x="-0.5" y="0.5"/>
>       <vertex index="8" x="-0.5" y="0.0"/>
>     </vertices>
>
>     <cells size="2">
>       <triangle index="0" affine="false" v0="0" v1="1" v2="2" v3="4"  
> v4="5" v5="6"/>
>       <triangle index="1" affine="true"  v0="3" v1="0" v2="2" v3="8"  
> v4="6" v5="7"/>
>     </cells>
>   </mesh>
> </dolfin>
>
> # note the variable `affine' that indicates whether the triangle is  
> actually higher order or not.
>
> ===============================================
>
> - Shawn
>
> On Tue, 19 Aug 2008, Anders Logg wrote:
>
>> On Tue, Aug 19, 2008 at 11:20:03AM -0400, Shawn Walker wrote:
>>> I think in dolfin's mesh.xml format, I will need to put the extra
>>> 'curved' vertices into a separate list of data.  If I put all of the (say
>>> 2nd order polynomial mesh) vertices into the usual list of vertices, then
>>> dolfin will think the matrix size is larger than it really is.  At least,
>>> that is what the output looks like.  Does this seem correct.
>>
>> One could either break the file format and introduce a <geometry> tag
>> that holds the geometry, or one could add a list of edge midpoints in
>> the data section. It could just be some <array> objects named "edge
>> midpoints x", "edge midpoints y" etc.
>>
>>> Also, I looked at this UFCCell.h file in Dolfin, and it looks like
>>> quadrilaterals and hexahedrons are not implemented.  Is this true?
>>
>> It's not supported anywhere except for in the UFC manual. :-)
>>
>>
>>
>>> - Shawn
>>>
>>> On Mon, 18 Aug 2008, Kent-Andre Mardal wrote:
>>>
>>>> On ma., 2008-08-18 at 09:01 -0400, Shawn Walker wrote:
>>>>> Sure, but I do want to have the option of some elements being curved and
>>>>> others not.  This is actually practical since only the elements near the
>>>>> boundary really need to be curved.
>>>>
>>>> Sure, I agree totally on this. But this requires more. Neither UFC, the
>>>> form compilers (FFC or SFC) nor the FEM basis function generators,  FIAT
>>>> and SyFi, support different polynomial degree at different places.
>>>>
>>>> Of course, in your case I guess it is only the geometrical mapping that
>>>> changes degree throughout the mesh and this may simplify alot eg. the
>>>> finite element basis functions will be the same throughout the mesh. The
>>>> dofs can also be reused.
>>>>
>>>> I only suggest starting with the implementation of a higher-order mesh
>>>> because it seems as a feasible first step.
>>>>
>>>> Kent
>>>>
>>>>
>>>>>
>>>>> Also, it is possible to evaluate some forms exactly, even if the local map
>>>>> is a polynomial.  Of course, this can be ignored for now...
>>>>>
>>>>> - Shawn
>>>>>
>>>>> On Mon, 18 Aug 2008, Kent-Andre Mardal wrote:
>>>>>
>>>>>> On s�008-08-17 at 00:08 -0400, Shawn Walker wrote:
>>>>>>> I think you are right.  I would propose the following:
>>>>>>>
>>>>>>> 1. Add a meshfunction to the mesh.xml file.  The name of it will
>>>>>>> correspond to a new enum variable in ufc::cell.  For example, say we have
>>>>>>> piecewise quadratics for the local map, so we call it "P2_map".  Maybe
>>>>>>> later this could be made a more permanent xml data type, but for now this
>>>>>>> should be ok.  The mesh function would just be a vector lagrange 2nd
>>>>>>> order function.  And there would also be a boolean function on each cell
>>>>>>> to say whether the element is straight or not.  This of course rules out
>>>>>>> mixing curved types, i.e. having P3 and P2 in the same mesh.  But I am ok
>>>>>>> with that for now.
>>>>>>>
>>>>>>> 2. In ufc::cell, we have another enum variable called "map_type":
>>>>>>>
>>>>>>> enum map_type {P1_map, P2_map, P3_map, etc...}
>>>>>>>
>>>>>>> Of course, this will only be for lagrange type polynomial maps.  But this
>>>>>>> enum variable can have other types (i.e. iso-geometric).  We then create a
>>>>>>> variable "map_type cell_map_type;"
>>>>>>
>>>>>> Yes, we probably would need one extra variable for the map.
>>>>>>
>>>>>> The way I see it, if one assume that all elements in the mesh has the
>>>>>> same order, the main work would be to create a dolfin mesh containing
>>>>>> a higher order geometry. Once this is done, you would need to ensure
>>>>>> that the coordinates in ufc::cell follows the UFC convention locally.
>>>>>> Finally, the computation of the Jacobian, the inverse of the Jacobian
>>>>>> and the determinant of the Jacobian has to be moved inside the
>>>>>> quadrature loop. This is the only modification to the form compiler, I
>>>>>> believe. The finite elements, its derivatives, dofs etc
>>>>>> stay the same.
>>>>>>
>>>>>> Hence, once a higher order mesh is up an running, there should only
>>>>>> be minimal modifications to the form compilers and UFC.
>>>>>>
>>>>>>
>>>>>>
>>>>>>>
>>>>>>> 3. When the mesh is read in, and the element is created, I think the
>>>>>>> "coordinates" variable can store the vertex positions of the element,
>>>>>>> including the extra ones for the curved sides.  If not, then some other
>>>>>>> variable would need to be created to store it; actually, that may be
>>>>>>> better.  Also, the "cell_map_type" will be set to either P1_map, or
>>>>>>> P2_map, or etc...
>>>>>>>
>>>>>>> 4. The tabulate_tensor routine will have a switch to pick the correct
>>>>>>> routine for computing the tensor coefs.  This will require modifying FFC
>>>>>>> which I don't know yet.
>>>>>>>
>>>>>>> 5. Another thing that should change will be the basis function
>>>>>>> evaluations.  This may be a pain.  Computing the inverse map to get the
>>>>>>> location of the coordinates in the reference element may be inconvenient.
>>>>>>> This can be ignored for now; this isn't really critical for what I would
>>>>>>> like to do.
>>>>>>>
>>>>>>> I am a little confused on where the ufc::cell data gets set by the data in
>>>>>>> the mesh.xml file.  I'm still very new to this code.  Any help would be
>>>>>>> appreciated.  If I am totally off on this, please say so!
>>>>>>>
>>>>>>
>>>>>> Have a look in UFCCell.
>>>>>>
>>>>>> Kent
>>>>>>
>>>>>>
>>>>>>
>>>>
>>>>
>>>>
>>
>>

Attachment: signature.asc
Description: Digital signature