diff --git a/Plugin/Gradient.cpp b/Plugin/Gradient.cpp deleted file mode 100644 index 422f62fb311c64cbc6774aac3a648cd6e03de00d..0000000000000000000000000000000000000000 --- a/Plugin/Gradient.cpp +++ /dev/null @@ -1,291 +0,0 @@ -// $Id: Gradient.cpp,v 1.3 2004-12-22 17:49:27 geuzaine Exp $ -// -// Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle -// -// This program is free software; you can redistribute it and/or modify -// it under the terms of the GNU General Public License as published by -// the Free Software Foundation; either version 2 of the License, or -// (at your option) any later version. -// -// This program is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -// GNU General Public License for more details. -// -// You should have received a copy of the GNU General Public License -// along with this program; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 -// USA. -// -// Please report all bugs and problems to <gmsh@geuz.org>. - -#include "Plugin.h" -#include "Gradient.h" -#include "List.h" -#include "Views.h" -#include "Context.h" -#include "Numeric.h" -#include <algorithm> -#include <functional> -#include <math.h> -#include <stdio.h> - -extern Context_T CTX; - -StringXNumber GradientOptions_Number[] = { - {GMSH_FULLRC, "iView", NULL, -1.} -}; - -extern "C" -{ - GMSH_Plugin *GMSH_RegisterGradientPlugin() - { - return new GMSH_GradientPlugin(); - } -} - -int GMSH_GradientPlugin::getNbOptions() const -{ - return sizeof(GradientOptions_Number) / sizeof(StringXNumber); -} - -StringXNumber *GMSH_GradientPlugin::getOption(int iopt) -{ - return &GradientOptions_Number[iopt]; -} - -GMSH_GradientPlugin::GMSH_GradientPlugin() -{ - ; -} - -void GMSH_GradientPlugin::getName(char *name) const -{ - strcpy(name, "Gradient"); -} - -void GMSH_GradientPlugin::getInfos(char *author, char *copyright, - char *help_text) const -{ - strcpy(author, "E. Marchandise"); - strcpy(copyright, "DGR (www.multiphysics.com)"); - strcpy(help_text, - "Plugin(Gradient) computes XXX\n" - "If `iView' < 0, the plugin is run\n" - "on the current view.\n" - "\n" - "Plugin(Gradient) creates one new view.\n"); -} - -void GMSH_GradientPlugin::catchErrorMessage(char *errorMessage) const -{ - strcpy(errorMessage, "Gradient failed..."); -} - -static void computeGradient(List_T *inList, int inNb, - List_T *outListVector, int *outNbVector, - int nbTime, int nbNod, int nbComp) -{ - if(!inNb) - return; - - int outNbComp, *outNb; - List_T *outList; - - outNbComp = 1; - outNb = outNbVector; - outList = outListVector; - - const int MAX_NOD = 4; - const int MAX_COMP= 3; - - int nb = List_Nbr(inList) / inNb; - - // on parcourt les elements - for(int i = 0; i < List_Nbr(inList); i += nb) { - - double *yy = (double *)List_Pointer_Fast(inList, i + nbNod); - - if(yy[0] > 0){ - for(int j = 0; j < 3 * nbNod; j++) - List_Add(outList, List_Pointer_Fast(inList, i + j)); - - for(int t = 0; t < nbTime; t++){ - - //on lit vecteurs positions et vitesses - //------------------------------------- - double *x = (double *)List_Pointer_Fast(inList, i); - double *y = (double *)List_Pointer_Fast(inList, i + nbNod); - double *z = (double *)List_Pointer_Fast(inList, i + 2 * nbNod); - double val[MAX_COMP][MAX_NOD]; - for(int k = 0; k < nbNod; k++){ - double *v = (double *)List_Pointer_Fast(inList, i + 3 * nbNod + - nbNod * nbComp * t + nbComp * k); - for(int l = 0; l < nbComp; l++){ - val[l][k] = v[l]; - } - } - - //on calcule le gradient des fonctions de forme - //--------------------------------------------- - double GradPhi_x[4][3]; - double GradPhi_ksi[4][3]; - double dx_dksi[3][3]; - double dksi_dx[3][3]; - double det; - - //triangles - //+++++++++ - if (nbNod == 3){ - double a[3], b[3], cross[3]; - a[0]= x[1]-x[0]; a[1]= y[1]-y[0]; a[2]= z[1]-z[0]; - b[0]= x[2]-x[0]; b[1]= y[2]-y[0]; b[2]= z[2]-z[0]; - prodve(a, b, cross); - dx_dksi[0][0] = x[1] - x[0]; dx_dksi[0][1] = x[2]-x[0]; dx_dksi[0][2] = cross[0]; - dx_dksi[1][0] = y[1] - y[0]; dx_dksi[1][1] = y[2]-y[0]; dx_dksi[1][2] = cross[1]; - dx_dksi[2][0] = z[1] - z[0]; dx_dksi[2][1] = z[2]-z[0]; dx_dksi[2][2] = cross[2]; - inv3x3(dx_dksi, dksi_dx, &det); - - GradPhi_ksi[0][0]= -1; GradPhi_ksi[0][1]= -1; GradPhi_ksi[0][2]= 0; - GradPhi_ksi[1][0]= 1; GradPhi_ksi[1][1]= 0; GradPhi_ksi[1][2]= 0; - GradPhi_ksi[2][0]= 0; GradPhi_ksi[2][1]= 1; GradPhi_ksi[2][2]= 0; - } - - //tetraedre - //++++++++++ - if (nbNod == 4){ - dx_dksi[0][0] = x[1] - x[0]; dx_dksi[0][1] = x[2]-x[0]; dx_dksi[0][2] = x[3]-x[0]; - dx_dksi[1][0] = y[1] - y[0]; dx_dksi[1][1] = y[2]-y[0]; dx_dksi[1][2] = y[3]-y[0]; - dx_dksi[2][0] = z[1] - z[0]; dx_dksi[2][1] = z[2]-z[0]; dx_dksi[2][2] = z[3]-z[0]; - inv3x3(dx_dksi, dksi_dx, &det); - GradPhi_ksi[0][0]= -1; GradPhi_ksi[0][1]= -1; GradPhi_ksi[0][2]= -1; - GradPhi_ksi[1][0]= 1; GradPhi_ksi[1][1]= 0; GradPhi_ksi[1][2]= 0; - GradPhi_ksi[2][0]= 0; GradPhi_ksi[2][1]= 1; GradPhi_ksi[2][2]= 0; - GradPhi_ksi[3][0]= 0; GradPhi_ksi[3][1]= 0; GradPhi_ksi[3][2]= 1; - } - - for(int k = 0; k < nbNod ; k++){ - for(int l = 0; l < 3; l++){ - GradPhi_x[k][l] = 0.0; - for(int m = 0; m < 3; m++){ - GradPhi_x[k][l] += GradPhi_ksi[k][m] * dksi_dx[l][m]; - } - } - } - - //on calcule le gradient des vitesses - //---------------------------------- - - //par noeud on a max 9 composantes - double res, grad[3][3], sym[3][3], asym[3][3]; - // for(int l = outNbComp; l < 9; l++) - // d[l] = 0.0; - for(int k = 0; k < nbComp; k++){ - for(int l = 0; l < 3; l++){ - grad[k][l] = 0.0; - for(int m = 0; m < nbNod; m++){ - grad[k][l] += val[k][m]* GradPhi_x[m][l]; - } - } - } - - for(int k=0; k<nbComp; k++){ - for(int l=0; l<3; l++){ - sym[k][l] = 0.5*(grad[k][l]+grad[l][k]); - asym[k][l] = 0.5*(grad[k][l]-grad[l][k]); - } - } - - double a[3][3]; - for(int i=0; i<3; i++) - for(int j=0; j<3; j++){ - a[i][j] = 0.0; - for(int k=0; k<3; k++) - a[i][j] += sym[i][k]*sym[k][j] + asym[i][k]*asym[k][j]; - } - - //calcul de lambda avec tri - double lambda[3]; - eigenvalue(a, lambda); - - for(int k = 0; k < nbNod; k++){ - for(int l = 0; l < outNbComp; l++){ - res = lambda[1]; - List_Add(outList, &res); - } - } - } - (*outNb)++; - } - } -} - -Post_View *GMSH_GradientPlugin::execute(Post_View * v) -{ - int iView = (int)GradientOptions_Number[0].def; - - if(iView < 0) - iView = v ? v->Index : 0; - - if(!List_Pointer_Test(CTX.post.list, iView)) { - Msg(GERROR, "View[%d] does not exist", iView); - return v; - } - - Post_View *v1 = *(Post_View **)List_Pointer(CTX.post.list, iView); - - Post_View *v2 = BeginView(1); - - // points - //-------- - // computeGradient( v1->SP, v1->NbSP, v2->VP, &v2->NbVP, v1->NbTimeStep, 1, 1); - // computeGradient( v1->VP, v1->NbVP, v2->TP, &v2->NbTP, v1->NbTimeStep, 1, 3); - - // lines - //------- - // computeGradient( v1->SL, v1->NbSL, v2->VL, &v2->NbVL, v1->NbTimeStep, 2, 1); - // computeGradient( v1->VL, v1->NbVL, v2->TL, &v2->NbTL, v1->NbTimeStep, 2, 3); - - // triangles - //------------ - // computeGradient( v1->ST, v1->NbST, v2->VT, &v2->NbVT, v1->NbTimeStep, 3, 1); - // computeGradient( v1->VT, v1->NbVT, v2->TT, &v2->NbTT, v1->NbTimeStep, 3, 3); - - // tets - //------ - // computeGradient( v1->SS, v1->NbSS, v2->VS, &v2->NbVS, v1->NbTimeStep, 4, 1); - computeGradient( v1->VS, v1->NbVS, v2->SS, &v2->NbSS, v1->NbTimeStep, 4, 3); - - /* - // quadrangles - //----------- - computeGradient( v1->SQ, v1->NbSQ, v2->VQ, &v2->NbVQ, v1->NbTimeStep, 4, 1); - computeGradient( v1->VQ, v1->NbVQ, v2->VQ, &v2->NbVQ, v1->NbTimeStep, 4, 3); - computeGradient( v1->TQ, v1->NbTQ, v2->VQ, &v2->NbVQ, v1->NbTimeStep, 4, 9); - // hexas - //-------- - computeGradient( v1->SH, v1->NbSH, v2->VH, &v2->NbVH, v1->NbTimeStep, 8, 1); - computeGradient( v1->VH, v1->NbVH, v2->VH, &v2->NbVH, v1->NbTimeStep, 8, 3); - computeGradient( v1->TH, v1->NbTH, v2->VH, &v2->NbVH, v1->NbTimeStep, 8, 9); - // prisms - //-------- - computeGradient( v1->SI, v1->NbSI, v2->VI, &v2->NbVI, v1->NbTimeStep, 6, 1); - computeGradient( v1->VI, v1->NbVI, v2->VI, &v2->NbVI, v1->NbTimeStep, 6, 3); - computeGradient( v1->TI, v1->NbTI, v2->VI, &v2->NbVI, v1->NbTimeStep, 6, 9); - // pyramids - //---------- - computeGradient( v1->SY, v1->NbSY, v2->VY, &v2->NbVY, v1->NbTimeStep, 5, 1); - computeGradient( v1->VY, v1->NbVY, v2->VY, &v2->NbVY, v1->NbTimeStep, 5, 3); - computeGradient( v1->TY, v1->NbTY, v2->VY, &v2->NbVY, v1->NbTimeStep, 5, 9); - */ - - // copy time data - for(int i = 0; i < List_Nbr(v1->Time); i++) - List_Add(v2->Time, List_Pointer(v1->Time, i)); - // finalize - char name[1024], filename[1024]; - sprintf(name, "%s_Gradient", v1->Name); - sprintf(filename, "%s_Gradient.pos", v1->Name); - EndView(v2, 1, filename, name); - return v2; -} diff --git a/Plugin/Gradient.h b/Plugin/Gradient.h deleted file mode 100644 index 516ef94946a69276f4dfe6ccecb859705fb8f163..0000000000000000000000000000000000000000 --- a/Plugin/Gradient.h +++ /dev/null @@ -1,43 +0,0 @@ -#ifndef _GRADIENT_H_ -#define _GRADIENT_H - -// Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle -// -// This program is free software; you can redistribute it and/or modify -// it under the terms of the GNU General Public License as published by -// the Free Software Foundation; either version 2 of the License, or -// (at your option) any later version. -// -// This program is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -// GNU General Public License for more details. -// -// You should have received a copy of the GNU General Public License -// along with this program; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 -// USA. -// -// Please report all bugs and problems to <gmsh@geuz.org>. - -#include "Plugin.h" -#include "Numeric.h" - -extern "C" -{ - GMSH_Plugin *GMSH_RegisterGradientPlugin(); -} - -class GMSH_GradientPlugin : public GMSH_Post_Plugin -{ -public: - GMSH_GradientPlugin(); - void getName(char *name) const; - void getInfos(char *author, char *copyright, char *helpText) const; - void catchErrorMessage(char *errorMessage) const; - int getNbOptions() const; - StringXNumber* getOption(int iopt); - Post_View *execute(Post_View *); -}; - -#endif diff --git a/Plugin/Lambda2.cpp b/Plugin/Lambda2.cpp index b968faad0fe9c5cdbcbf0ee19f1f2f857557c54e..7cd6bfe22219defc40982ba84da09c97a4cd26e5 100644 --- a/Plugin/Lambda2.cpp +++ b/Plugin/Lambda2.cpp @@ -1,4 +1,4 @@ -// $Id: Lambda2.cpp,v 1.3 2004-12-07 04:52:27 geuzaine Exp $ +// $Id: Lambda2.cpp,v 1.4 2004-12-27 19:18:12 geuzaine Exp $ // // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle // @@ -25,8 +25,7 @@ #include "Views.h" #include "Context.h" #include "Numeric.h" -#include <algorithm> -#include <functional> +#include "ShapeFunctions.h" #include <math.h> #include <stdio.h> @@ -61,19 +60,22 @@ void GMSH_Lambda2Plugin::getInfos(char *author, char *copyright, strcpy(author, "E. Marchandise"); strcpy(copyright, "DGR (www.multiphysics.com)"); strcpy(help_text, - "Plugin(Lambda2) computes\n" - "the eigenvalues Lambda (1,2,3)\n" - "of the tensor S_ik S_kj + Om_ik Om_kj,\n" - "where S_ij = 0.5 (ui,j + uj,i)) and\n" - "Om_ij = 0.5 (ui,j - uj,i) are respectively\n" - "the symmetric and antisymmetric parts\n" - "of the velocity gradient tensor.\n" - "Vortices are well-represented by regions\n" - "where Lambda (2) is negative.\n" - "If `iView' < 0, the plugin is run\n" - "on the current view.\n" + "Plugin(Lambda2) computes the eigenvalues\n" + "Lambda(1,2,3) of the tensor (S_ik S_kj +\n" + "Om_ik Om_kj), where S_ij = 0.5 (ui,j + uj,i)\n" + "and Om_ij = 0.5 (ui,j - uj,i) are respectively\n" + "the symmetric and antisymmetric parts of the\n" + "velocity gradient tensor. Vortices are well-\n" + "represented by regions where Lambda(2) is\n" + "negative. If `iView' contains tensor elements,\n" + "the plugin directly uses the tensors as the\n" + "values of the velocity gradient tensor; if\n" + "`iView' contains vector elements, the plugin\n" + "uses them as the velocities from which to derive\n" + "the velocity gradient tensor. If `iView' < 0,\n" + "the plugin is run on the current view.\n" "\n" - "Plugin(Lambda2) creates one noew view.\n"); + "Plugin(Lambda2) creates one new view.\n"); } int GMSH_Lambda2Plugin::getNbOptions() const @@ -92,69 +94,131 @@ void GMSH_Lambda2Plugin::catchErrorMessage(char *errorMessage) const } static void eigen(List_T *inList, int inNb, - List_T *outListScalar, int *outNbScalar, + List_T *outList, int *outNb, int nbTime, int nbNod, int nbComp, int lam) { - if(!inNb) + if(!inNb || (nbComp != 3 && nbComp != 9) || lam < 1 || lam > 3) return; - - int outNbComp, *outNb; - List_T *outList; - - outNbComp = 1; - outNb = outNbScalar; - outList = outListScalar; - + + // loop on elements int nb = List_Nbr(inList) / inNb; - - //boucle sur les elements - //----------------------- for(int i = 0; i < List_Nbr(inList); i += nb) { - + + // copy node coordinates for(int j = 0; j < 3 * nbNod; j++) List_Add(outList, List_Pointer_Fast(inList, i + j)); + // loop on time steps for(int j = 0; j < nbTime; j++){ - - double *v = (double *)List_Pointer_Fast(inList, i + 3 * nbNod + - nbNod * nbComp * j + nbComp * 0); + + double *x = (double *)List_Pointer_Fast(inList, i); + double *y = (double *)List_Pointer_Fast(inList, i + nbNod); + double *z = (double *)List_Pointer_Fast(inList, i + 2 * nbNod); + double GradVel[3][3]; - GradVel[0][0] = v[0]; GradVel[0][1] = v[1]; GradVel[0][2] = v[2]; - GradVel[1][0] = v[3]; GradVel[1][1] = v[4]; GradVel[1][2] = v[5]; - GradVel[2][0] = v[6]; GradVel[2][1] = v[7]; GradVel[2][2] = v[8]; - - //par element - //------------ - - //calcul partie sym et antisymetrique + + if(nbComp == 9){ + // val is the velocity gradient tensor, we assume that is is + // contant per element + double *v = (double *)List_Pointer_Fast(inList, i + 3 * nbNod + + nbNod * nbComp * j + nbComp * 0); + GradVel[0][0] = v[0]; GradVel[0][1] = v[1]; GradVel[0][2] = v[2]; + GradVel[1][0] = v[3]; GradVel[1][1] = v[4]; GradVel[1][2] = v[5]; + GradVel[2][0] = v[6]; GradVel[2][1] = v[7]; GradVel[2][2] = v[8]; + } + else if(nbComp == 3){ + // FIXME: the following could be greatly simplified and + // generalized by using the classes in ShapeFunctions.h! + + // val contains the velocities: compute the gradient tensor + // from them + const int MAX_NOD = 4; + const int MAX_COMP= 3; + double val[MAX_COMP][MAX_NOD]; + for(int k = 0; k < nbNod; k++){ + double *v = (double *)List_Pointer_Fast(inList, i + 3 * nbNod + + nbNod * nbComp * j + nbComp * k); + for(int l = 0; l < 3; l++){ + val[l][k] = v[l]; + } + } + // compute gradient of shape functions + double GradPhi_x[4][3]; + double GradPhi_ksi[4][3]; + double dx_dksi[3][3]; + double dksi_dx[3][3]; + double det; + if(nbNod == 3){ // triangles + double a[3], b[3], cross[3]; + a[0]= x[1]-x[0]; a[1]= y[1]-y[0]; a[2]= z[1]-z[0]; + b[0]= x[2]-x[0]; b[1]= y[2]-y[0]; b[2]= z[2]-z[0]; + prodve(a, b, cross); + dx_dksi[0][0] = x[1] - x[0]; dx_dksi[0][1] = x[2]-x[0]; dx_dksi[0][2] = cross[0]; + dx_dksi[1][0] = y[1] - y[0]; dx_dksi[1][1] = y[2]-y[0]; dx_dksi[1][2] = cross[1]; + dx_dksi[2][0] = z[1] - z[0]; dx_dksi[2][1] = z[2]-z[0]; dx_dksi[2][2] = cross[2]; + inv3x3(dx_dksi, dksi_dx, &det); + GradPhi_ksi[0][0]= -1; GradPhi_ksi[0][1]= -1; GradPhi_ksi[0][2]= 0; + GradPhi_ksi[1][0]= 1; GradPhi_ksi[1][1]= 0; GradPhi_ksi[1][2]= 0; + GradPhi_ksi[2][0]= 0; GradPhi_ksi[2][1]= 1; GradPhi_ksi[2][2]= 0; + } + else if (nbNod == 4){ // tetraedre + dx_dksi[0][0] = x[1] - x[0]; dx_dksi[0][1] = x[2]-x[0]; dx_dksi[0][2] = x[3]-x[0]; + dx_dksi[1][0] = y[1] - y[0]; dx_dksi[1][1] = y[2]-y[0]; dx_dksi[1][2] = y[3]-y[0]; + dx_dksi[2][0] = z[1] - z[0]; dx_dksi[2][1] = z[2]-z[0]; dx_dksi[2][2] = z[3]-z[0]; + inv3x3(dx_dksi, dksi_dx, &det); + GradPhi_ksi[0][0]= -1; GradPhi_ksi[0][1]= -1; GradPhi_ksi[0][2]= -1; + GradPhi_ksi[1][0]= 1; GradPhi_ksi[1][1]= 0; GradPhi_ksi[1][2]= 0; + GradPhi_ksi[2][0]= 0; GradPhi_ksi[2][1]= 1; GradPhi_ksi[2][2]= 0; + GradPhi_ksi[3][0]= 0; GradPhi_ksi[3][1]= 0; GradPhi_ksi[3][2]= 1; + } + else{ + Msg(GERROR, "Lambda2 not ready for this type of element"); + return; + } + for(int k = 0; k < nbNod ; k++){ + for(int l = 0; l < 3; l++){ + GradPhi_x[k][l] = 0.0; + for(int m = 0; m < 3; m++){ + GradPhi_x[k][l] += GradPhi_ksi[k][m] * dksi_dx[l][m]; + } + } + } + // compute gradient of velocities + for(int k = 0; k < 3; k++){ + for(int l = 0; l < 3; l++){ + GradVel[k][l] = 0.0; + for(int m = 0; m < nbNod; m++){ + GradVel[k][l] += val[k][m]* GradPhi_x[m][l]; + } + } + } + } + + // compute the sym and antisymetric parts double sym[3][3]; double asym[3][3]; - for(int i=0; i<3; i++){ - for(int j=0; j<3; j++){ - sym[i][j] = 0.5*(GradVel[i][j]+GradVel[j][i]); - asym[i][j] = 0.5*(GradVel[i][j]-GradVel[j][i]); - } + for(int m = 0; m < 3; m++){ + for(int n = 0; n < 3; n++){ + sym[m][n] = 0.5 * (GradVel[m][n] + GradVel[n][m]); + asym[m][n] = 0.5 * (GradVel[m][n] - GradVel[n][m]); + } } double a[3][3]; - for(int i=0; i<3; i++) - for(int j=0; j<3; j++){ - a[i][j] = 0.0; - for(int k=0; k<3; k++) - a[i][j] += sym[i][k]*sym[k][j] + asym[i][k]*asym[k][j]; - } + for(int m = 0; m < 3; m++){ + for(int n = 0; n < 3; n++){ + a[m][n] = 0.0; + for(int l = 0; l < 3; l++) + a[m][n] += sym[m][l] * sym[l][n] + asym[m][l] * asym[l][n]; + } + } - //calcul de lambda avec tri + // compute the eigenvalues double lambda[3]; eigenvalue(a, lambda); - - double res; - for(int k = 0; k < nbNod; k++){ - for(int l = 0; l < outNbComp; l++){ - res = lambda[lam-1]; - List_Add(outList, &res); - } - } + for(int k = 0; k < nbNod; k++) + List_Add(outList, &lambda[lam-1]); } + (*outNb)++; } } @@ -174,17 +238,35 @@ Post_View *GMSH_Lambda2Plugin::execute(Post_View * v) Post_View *v1 = *(Post_View **)List_Pointer(CTX.post.list, iView); Post_View *v2 = BeginView(1); - - eigen(v1->TT, v1->NbTT, v2->ST, &v2->NbSS, v1->NbTimeStep, 3, 9, ev); - eigen(v1->TS, v1->NbTS, v2->SS, &v2->NbSS, v1->NbTimeStep, 4, 9, ev); + // assume that the tensors contain the velocity gradient tensor + eigen(v1->TP, v1->NbTP, v2->SP, &v2->NbSP, v1->NbTimeStep, 1, 9, ev); + eigen(v1->TL, v1->NbTL, v2->SL, &v2->NbSL, v1->NbTimeStep, 2, 9, ev); + eigen(v1->TT, v1->NbTT, v2->ST, &v2->NbST, v1->NbTimeStep, 3, 9, ev); + eigen(v1->TQ, v1->NbTQ, v2->SQ, &v2->NbSQ, v1->NbTimeStep, 4, 9, ev); + eigen(v1->TS, v1->NbTS, v2->SS, &v2->NbSS, v1->NbTimeStep, 4, 9, ev); + eigen(v1->TH, v1->NbTH, v2->SH, &v2->NbSH, v1->NbTimeStep, 8, 9, ev); + eigen(v1->TI, v1->NbTI, v2->SI, &v2->NbSI, v1->NbTimeStep, 6, 9, ev); + eigen(v1->TY, v1->NbTY, v2->SY, &v2->NbSY, v1->NbTimeStep, 5, 9, ev); + + // assume that the vectors contain the velocities + // FIXME: only implemented for tri/tet at the moment + //eigen(v1->VP, v1->NbVP, v2->SP, &v2->NbSP, v1->NbTimeStep, 1, 3, ev); + //eigen(v1->VL, v1->NbVL, v2->SL, &v2->NbSL, v1->NbTimeStep, 2, 3, ev); + eigen(v1->VT, v1->NbVT, v2->ST, &v2->NbST, v1->NbTimeStep, 3, 3, ev); + //eigen(v1->VQ, v1->NbVQ, v2->SQ, &v2->NbSQ, v1->NbTimeStep, 4, 3, ev); + eigen(v1->VS, v1->NbVS, v2->SS, &v2->NbSS, v1->NbTimeStep, 4, 3, ev); + //eigen(v1->VH, v1->NbVH, v2->SH, &v2->NbSH, v1->NbTimeStep, 8, 3, ev); + //eigen(v1->VI, v1->NbVI, v2->SI, &v2->NbSI, v1->NbTimeStep, 6, 3, ev); + //eigen(v1->VY, v1->NbVY, v2->SY, &v2->NbSY, v1->NbTimeStep, 5, 3, ev); + // copy time data for(int i = 0; i < List_Nbr(v1->Time); i++) List_Add(v2->Time, List_Pointer(v1->Time, i)); // finalize char name[1024], filename[1024]; - sprintf(name, "%s_Extract", v1->Name); - sprintf(filename, "%s_Extract.pos", v1->Name); + sprintf(name, "%s_Lambda2", v1->Name); + sprintf(filename, "%s_Lambda2.pos", v1->Name); EndView(v2, 1, filename, name); return v2; } diff --git a/Plugin/Lambda2.h b/Plugin/Lambda2.h index f7d9f931744cfed172980185e8a7762a6f3a7a00..1b391f9ccd85e5f37e10fd10d79ba79f5a7286e7 100644 --- a/Plugin/Lambda2.h +++ b/Plugin/Lambda2.h @@ -44,9 +44,6 @@ public: StringXNumber* getOption(int iopt); Post_View *execute(Post_View *); - - - -//--------------------------------------------------------- }; + #endif diff --git a/Plugin/Makefile b/Plugin/Makefile index 36946e8e334c2566613e4a1df3de3eeb8fad0820..204a348613f490c62b515e13b1c355abb84c9074 100644 --- a/Plugin/Makefile +++ b/Plugin/Makefile @@ -1,4 +1,4 @@ -# $Id: Makefile,v 1.69 2004-12-27 09:17:44 geuzaine Exp $ +# $Id: Makefile,v 1.70 2004-12-27 19:18:12 geuzaine Exp $ # # Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle # @@ -30,7 +30,7 @@ SRC = Plugin.cpp\ Levelset.cpp\ CutPlane.cpp CutSphere.cpp CutMap.cpp \ Smooth.cpp CutParametric.cpp\ - Gradient.cpp Lambda2.cpp\ + Lambda2.cpp\ Eigenvectors.cpp\ Octree.cpp OctreeInternals.cpp OctreePost.cpp\ StreamLines.cpp CutGrid.cpp\ @@ -82,8 +82,8 @@ Plugin.o: Plugin.cpp Plugin.h ../Common/Options.h ../Common/Message.h \ StructuralSolver.h ../Geo/Geo.h ../Mesh/Mesh.h ../Mesh/Vertex.h \ ../Mesh/Element.h ../Mesh/Simplex.h ../Mesh/Face.h ../Mesh/Edge.h \ ../Geo/ExtrudeParams.h ../Mesh/STL.h ../Mesh/Metric.h ../Mesh/Matrix.h \ - ../Common/GmshUI.h Eigenvectors.h Evaluate.h Probe.h \ - ../Common/Context.h + ../Common/GmshUI.h Eigenvectors.h Lambda2.h ../Numeric/Numeric.h \ + Evaluate.h Probe.h ../Common/Context.h Levelset.o: Levelset.cpp Levelset.h Plugin.h ../Common/Options.h \ ../Common/Message.h ../Common/Views.h ../Common/ColorTable.h \ ../DataStr/List.h ../Common/VertexArray.h ../Common/SmoothNormals.h \ @@ -122,16 +122,11 @@ CutParametric.o: CutParametric.cpp OctreePost.h Octree.h \ ../Common/Message.h ../Common/Views.h ../Common/ColorTable.h \ ../DataStr/List.h ../Common/VertexArray.h ../Common/SmoothNormals.h \ ../Common/GmshMatrix.h ../Common/AdaptiveViews.h ../Common/Context.h -Gradient.o: Gradient.cpp Plugin.h ../Common/Options.h ../Common/Message.h \ - ../Common/Views.h ../Common/ColorTable.h ../DataStr/List.h \ - ../Common/VertexArray.h ../Common/SmoothNormals.h \ - ../Common/GmshMatrix.h ../Common/AdaptiveViews.h Gradient.h \ - ../Numeric/Numeric.h ../Common/Context.h Lambda2.o: Lambda2.cpp Plugin.h ../Common/Options.h ../Common/Message.h \ ../Common/Views.h ../Common/ColorTable.h ../DataStr/List.h \ ../Common/VertexArray.h ../Common/SmoothNormals.h \ ../Common/GmshMatrix.h ../Common/AdaptiveViews.h Lambda2.h \ - ../Numeric/Numeric.h ../Common/Context.h + ../Numeric/Numeric.h ../Common/Context.h ShapeFunctions.h Eigenvectors.o: Eigenvectors.cpp Plugin.h ../Common/Options.h \ ../Common/Message.h ../Common/Views.h ../Common/ColorTable.h \ ../DataStr/List.h ../Common/VertexArray.h ../Common/SmoothNormals.h \ diff --git a/Plugin/Plugin.cpp b/Plugin/Plugin.cpp index e59761ce81abed42a910e0a29e49a59fa88a7e07..ecdb10e2aa1ca3c74a7c3708a61165f2cebf050e 100644 --- a/Plugin/Plugin.cpp +++ b/Plugin/Plugin.cpp @@ -1,4 +1,4 @@ -// $Id: Plugin.cpp,v 1.69 2004-12-27 07:31:14 geuzaine Exp $ +// $Id: Plugin.cpp,v 1.70 2004-12-27 19:18:12 geuzaine Exp $ // // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle // @@ -52,6 +52,7 @@ #include "DisplacementRaise.h" #include "StructuralSolver.h" #include "Eigenvectors.h" +#include "Lambda2.h" #include "Evaluate.h" #include "Probe.h" #include "Context.h" @@ -204,6 +205,8 @@ void GMSH_PluginManager::registerDefaultPlugins() ("Remove", GMSH_RegisterRemovePlugin())); allPlugins.insert(std::pair < char *, GMSH_Plugin * > ("Eigenvectors", GMSH_RegisterEigenvectorsPlugin())); + allPlugins.insert(std::pair < char *, GMSH_Plugin * > + ("Lambda2", GMSH_RegisterLambda2Plugin())); allPlugins.insert(std::pair < char *, GMSH_Plugin * > ("Probe", GMSH_RegisterProbePlugin())); #if defined(HAVE_TRIANGLE) diff --git a/doc/texinfo/opt_plugin.texi b/doc/texinfo/opt_plugin.texi index e2e1762ce7ee0b79111c7c7f7df9949011bab0e9..877cc7c7f3e743afc73c06307aefcf5e91e33f45 100644 --- a/doc/texinfo/opt_plugin.texi +++ b/doc/texinfo/opt_plugin.texi @@ -352,6 +352,32 @@ Default value: @code{-1} Default value: @code{0} @end table +@item Plugin(Lambda2) +Plugin(Lambda2) computes the eigenvalues +Lambda(1,2,3) of the tensor (S_ik S_kj + +Om_ik Om_kj), where S_ij = 0.5 (ui,j + uj,i) +and Om_ij = 0.5 (ui,j - uj,i) are respectively +the symmetric and antisymmetric parts of the +velocity gradient tensor. Vortices are well- +represented by regions where Lambda(2) is +negative. If `iView' contains tensor elements, +the plugin directly uses the tensors as the +values of the velocity gradient tensor; if +`iView' contains vector elements, the plugin +uses them as the velocities from which to derive +the velocity gradient tensor. If `iView' < 0, +the plugin is run on the current view. + +Plugin(Lambda2) creates one new view. + +Numeric options: +@table @code +@item Eigenvalue +Default value: @code{2} +@item iView +Default value: @code{-1} +@end table + @item Plugin(Probe) Plugin(Probe) gets the value of the simplectic view `iView' at the point (`X',`Y',`Z'). If `iView' < 0,