Subject: [FFC-dev] Re: [DOLFIN-dev] Updated Tecplot interface

From: Anders Logg <[email-masked]>
Date: Tue, 5 Jul 2005 12:20:01 -0500

On Tue, Jul 05, 2005 at 06:40:15PM +0200, Garth N. Wells wrote:
> Attached are updated Tecplot files (src/kernel/io/TecplotFile.cpp and
> src/kernel/io/TecplotFile.h). 

Thanks. I will add before the new release.

> I'm trying to stop with Tecplot for postprocessing and have started using Paraview
> (www.paraview.org) which is open source. If there is interest, I'd be happy to
> contribute the DOLFIN interface.

I would definitely be interested in having a DOLFIN interface for Paraview.

> Is it possible to compute second derivatives of the basis functions with FFC +
> FIAT? Now that higher-order elements can be implemented, I'd find this very useful.
>
> Garth

Yes, you can compute second derivatives (or any other derivatives),
but there is no good way yet to postprocess a solution and compute
it's second derivatives.

What you can do, is to define a linear system that computes the
projection of the second derivatives on to some space. As an example,
the following defines a linear system for computing the projection of
the x-derivative of a piecewise quadratic function u onto piecewise
linears:

P1 = FiniteElement("Lagrange", "triangle", 1)
P2 = FiniteElement("Lagrange", "triangle", 2)

v  = BasisFunction(P1)
ux = BasisFunction(P1)
u  = Function(P2)

a = v * ux * dx
L = v * u.dx(0) * dx

/Anders

> // Copyright (C) 2004 Harald Svensson.
> // Licensed under the GNU GPL Version 2.
> //
> // Modified by Anders Logg, 2004.
> // Modified by Garth N. Wells, 2005.
> 
> #include <dolfin/Mesh.h>
> #include <dolfin/TecplotFile.h>
> #include <dolfin/FiniteElement.h>
> 
> 
> using namespace dolfin;
> 
> //-?---------------------------------------------------------------------------
> TecplotFile::TecplotFile(const std::string filename) : GenericFile(filename)
> {
>   type = "TECPLOT";
> }
> //-?---------------------------------------------------------------------------
> TecplotFile::~TecplotFile()
> {
>   // Do nothing
> }
> //-?---------------------------------------------------------------------------
> void TecplotFile::operator<<(Mesh& mesh)
> {
> 
>   dolfin_info("Saving mesh to Tecplot file.");
> 
>   // Open file
>   FILE* fp = fopen(filename.c_str(), "a");
> 
> 	// Write header
>   fprintf(fp, "TITLE = \"Dolfin output\"  \n");
>   fprintf(fp, "VARIABLES = ");
>   if ( mesh.type() == Mesh::tetrahedra ){
> 	  fprintf(fp, " X1  X2  X3 \n");
> 	  fprintf(fp, "ZONE T = \" - \" N = %8d, E = %8d, DATAPACKING = POINT, ZONETYPE=FETETRAHEDRON \n", mesh.noNodes(), mesh.noCells());
>   }
> 	if ( mesh.type() == Mesh::triangles ){    
> 	  fprintf(fp, " X1  X2  X3 \n");
>     fprintf(fp, "ZONE T = \" - \"  N = %8d, E = %8d,  DATAPACKING = POINT, ZONETYPE=FETRIANGLE \n",   mesh.noNodes(), mesh.noCells());
>    }
> 
>   // Write node locations
>   for (NodeIterator n(mesh); !n.end(); ++n)
>   {
>     Point   p = n->coord();
> 
>     if ( mesh.type() == Mesh::tetrahedra )  fprintf(fp," %e %e %e \n",p.x, p.y, p.z);
>     if ( mesh.type() == Mesh::triangles )     fprintf(fp," %e %e  ",p.x, p.y);
>     fprintf(fp,"\n");
> 
>   }
> 
>   // Write cell connectivity
>   for (CellIterator c(mesh); !c.end(); ++c)
>   {
>     for (NodeIterator n(c); !n.end(); ++n) fprintf(fp," %8d ",n->id()+1);
>     fprintf(fp," \n");
>   }  
> 
>   // Close file
>   fclose(fp);
> 
> }
> //-?---------------------------------------------------------------------------
> void TecplotFile::operator<<(Function& u)
> {
> 
>   const FiniteElement& element = u.element();
> 	
> 	FILE *fp = fopen(filename.c_str(), "a");
>     
>   uint ShapeDim  = element.shapedim();
>   uint VectorDim = 0;
> 
>   if ( element.rank() == 0 )
>   {
>     VectorDim = 1;
>   }
>   else if ( element.rank() == 1 )
>   {
>     VectorDim = element.tensordim(0);
>   }
> 	
>   // Write mesh the first time
>   if ( u.number() == 0 )
>   {
> 	  // Write header
>     fprintf(fp, "TITLE = \"Dolfin output\"  \n");
>     fprintf(fp, "VARIABLES = ");
>     for (uint i=0; i<ShapeDim; ++i)   fprintf(fp, " X%d  ", i+1);	  
>     for (uint i=0; i<VectorDim; ++i)  fprintf(fp, " U%d  ", i+1);	  
>     fprintf(fp, "\n");	  
>     if ( u.mesh().type() == Mesh::tetrahedra )
> 	     fprintf(fp, "ZONE T = \"%6d\" N = %8d, E = %8d, DATAPACKING = POINT, ZONETYPE=FETETRAHEDRON \n", u.number()+1, u.mesh().noNodes(), u.mesh().noCells());
>     if ( u.mesh().type() == Mesh::triangles )
>        fprintf(fp, "ZONE T = \"%6d\"  N = %8d, E = %8d, DATAPACKING = POINT, ZONETYPE=FETRIANGLE \n",   u.number()+1, u.mesh().noNodes(), u.mesh().noCells());
> 
> 
>     // Write node locations and results
>     for (NodeIterator n(u.mesh()); !n.end(); ++n)
>     {
>       Point p = n->coord();
> 
>       if ( u.mesh().type() == Mesh::tetrahedra )  fprintf(fp," %e %e %e \n", p.x, p.y, p.z);
>       if ( u.mesh().type() == Mesh::triangles )     fprintf(fp," %e %e  ", p.x, p.y);
>       for (uint i=0; i < VectorDim; ++i) fprintf(fp,"%e ", u(*n,i) );
>       fprintf(fp,"\n");
> 
>       }
> 
>       // Write cell connectivity
>      for (CellIterator c(u.mesh()); !c.end(); ++c)
>      {
>        for (NodeIterator n(c); !n.end(); ++n) fprintf(fp," %8d ",n->id()+1);
>        fprintf(fp," \n");
>      }  
> 
>    }
> 	  
> 
>   // Write data for seccond and subsequent times
>   if ( u.number() != 0 )
>   {
>       // Write header
> 	  if ( u.mesh().type() == Mesh::tetrahedra )
>       	   fprintf(fp, "ZONE T = \"%6d\" N = %8d, E = %8d,  DATAPACKING = POINT, ZONETYPE=FETETRAHEDRON, VARSHARELIST = ([1-3]=1) CONNECTIVITYSHAREZONE=1 \n", u.number()+1, u.mesh().noNodes(), u.mesh().noCells());
>      if ( u.mesh().type() == Mesh::triangles )
>            fprintf(fp, "ZONE T = \"%6d\"  N = %8d, E = %8d, DATAPACKING = POINT, ZONETYPE=FETRIANGLE, VARSHARELIST = ([1,2]=1) CONNECTIVITYSHAREZONE=1 \n", u.number()+1, u.mesh().noNodes(), u.mesh().noCells());
> 
> 
>       // Write node locations and results
>       for (NodeIterator n(u.mesh()); !n.end(); ++n)
>       {
> 
>         for (uint i=0; i < VectorDim; ++i) fprintf(fp,"%e ", u(*n,i) );
>         fprintf(fp,"\n");
> 
>       }
>   }
>     
>   // Close file
>   fclose(fp);
>   
>   // Increase the number of times we have saved the function
>   ++u;
> 
>   cout << "Saved function " << u.name() << " (" << u.label()
>        << ") to file " << filename << " in Tecplot format." << endl;
> 
> 
> 
> }
> //-?---------------------------------------------------------------------------

> // Copyright (C) 2004 Harald Svensson.
> // Licensed under the GNU GPL Version 2.
> //
> // Modified by Garth N. Wells, 2005.
> // Modified by Anders Logg, 2005.
> // Modified by Garth N. Wells, 2005.
> 
> #ifndef __TECPLOT_FILE_H
> #define __TECPLOT_FILE_H
> 
> #include <dolfin/GenericFile.h>
> 
> namespace dolfin
> {
> 
>   class TecplotFile : public GenericFile
>   {
>   public:
> 
>     TecplotFile(const std::string filename);
>     ~TecplotFile();
> 
>     // Input
> 
>     // Output
>     void operator<< (Mesh& mesh);
>     void operator<< (Function& u);
> 
>   };
>   
> }
> 
> #endif

> _______________________________________________
> DOLFIN-dev mailing list
> [email-masked]
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev

-- 
Anders Logg
Research Assistant Professor
Toyota Technological Institute at Chicago
http://www.tti-c.org/logg/