1 // Gmsh - Copyright (C) 1997-2021 C. Geuzaine, J.-F. Remacle
2 //
3 // See the LICENSE.txt file in the Gmsh root directory for license information.
4 // Please report all issues on https://gitlab.onelab.info/gmsh/gmsh/issues.
5
6 #include "Eigenvalues.h"
7 #include "Numeric.h"
8 #include "GmshDefines.h"
9
10 StringXNumber EigenvaluesOptions_Number[] = {
11 {GMSH_FULLRC, "View", nullptr, -1.}};
12
13 extern "C" {
GMSH_RegisterEigenvaluesPlugin()14 GMSH_Plugin *GMSH_RegisterEigenvaluesPlugin()
15 {
16 return new GMSH_EigenvaluesPlugin();
17 }
18 }
19
getHelp() const20 std::string GMSH_EigenvaluesPlugin::getHelp() const
21 {
22 return "Plugin(Eigenvalues) computes the three real "
23 "eigenvalues of each tensor in the view `View'.\n\n"
24 "If `View' < 0, the plugin is run on the current view.\n\n"
25 "Plugin(Eigenvalues) creates three new list-based scalar views.";
26 }
27
getNbOptions() const28 int GMSH_EigenvaluesPlugin::getNbOptions() const
29 {
30 return sizeof(EigenvaluesOptions_Number) / sizeof(StringXNumber);
31 }
32
getOption(int iopt)33 StringXNumber *GMSH_EigenvaluesPlugin::getOption(int iopt)
34 {
35 return &EigenvaluesOptions_Number[iopt];
36 }
37
execute(PView * v)38 PView *GMSH_EigenvaluesPlugin::execute(PView *v)
39 {
40 int iView = (int)EigenvaluesOptions_Number[0].def;
41
42 PView *v1 = getView(iView, v);
43 if(!v1) return v;
44
45 PViewData *data1 = getPossiblyAdaptiveData(v1);
46 if(data1->hasMultipleMeshes()) {
47 Msg::Error("Eigenvalues plugin cannot be run on multi-mesh views");
48 return v;
49 }
50
51 PView *min = new PView();
52 PView *mid = new PView();
53 PView *max = new PView();
54
55 PViewDataList *dmin = getDataList(min);
56 PViewDataList *dmid = getDataList(mid);
57 PViewDataList *dmax = getDataList(max);
58
59 for(int ent = 0; ent < data1->getNumEntities(0); ent++) {
60 for(int ele = 0; ele < data1->getNumElements(0, ent); ele++) {
61 if(data1->skipElement(0, ent, ele)) continue;
62 int numComp = data1->getNumComponents(0, ent, ele);
63 if(numComp != 9) continue;
64 int type = data1->getType(0, ent, ele);
65 int numNodes = data1->getNumNodes(0, ent, ele);
66 std::vector<double> *outmin = dmin->incrementList(1, type, numNodes);
67 std::vector<double> *outmid = dmid->incrementList(1, type, numNodes);
68 std::vector<double> *outmax = dmax->incrementList(1, type, numNodes);
69 if(!outmin || !outmid || !outmax) continue;
70 double xyz[3][8];
71 for(int nod = 0; nod < numNodes; nod++)
72 data1->getNode(0, ent, ele, nod, xyz[0][nod], xyz[1][nod], xyz[2][nod]);
73 for(int i = 0; i < 3; i++) {
74 for(int nod = 0; nod < numNodes; nod++) {
75 outmin->push_back(xyz[i][nod]);
76 outmid->push_back(xyz[i][nod]);
77 outmax->push_back(xyz[i][nod]);
78 }
79 }
80 for(int step = 0; step < data1->getNumTimeSteps(); step++) {
81 for(int nod = 0; nod < numNodes; nod++) {
82 double val[9], w[3];
83 for(int comp = 0; comp < numComp; comp++)
84 data1->getValue(step, ent, ele, nod, comp, val[comp]);
85 double A[3][3] = {{val[0], val[1], val[2]},
86 {val[3], val[4], val[5]},
87 {val[6], val[7], val[8]}};
88 eigenvalue(A, w);
89 outmin->push_back(w[2]);
90 outmid->push_back(w[1]);
91 outmax->push_back(w[0]);
92 }
93 }
94 }
95 }
96
97 for(int i = 0; i < data1->getNumTimeSteps(); i++) {
98 double time = data1->getTime(i);
99 dmin->Time.push_back(time);
100 dmid->Time.push_back(time);
101 dmax->Time.push_back(time);
102 }
103 dmin->setName(data1->getName() + "_MinEigenvalues");
104 dmin->setFileName(data1->getName() + "_MinEigenvalues.pos");
105 dmin->finalize();
106 dmid->setName(data1->getName() + "_MidEigenvalues");
107 dmid->setFileName(data1->getName() + "_MidEigenvalues.pos");
108 dmid->finalize();
109 dmax->setName(data1->getName() + "_MaxEigenvalues");
110 dmax->setFileName(data1->getName() + "_MaxEigenvalues.pos");
111 dmax->finalize();
112
113 return nullptr;
114 }
115