From 934bf0c1c06a1a02701ebddc5e17e7a5ae7b19e9 Mon Sep 17 00:00:00 2001 From: Christophe Geuzaine <cgeuzaine@ulg.ac.be> Date: Wed, 31 Oct 2012 09:05:37 +0000 Subject: [PATCH] pp --- Numeric/GaussQuadraturePri.cpp | 44 ++++++++++++++++------------------ Numeric/fullMatrix.h | 15 ++++++++---- doc/gmsh.html | 3 ++- 3 files changed, 33 insertions(+), 29 deletions(-) diff --git a/Numeric/GaussQuadraturePri.cpp b/Numeric/GaussQuadraturePri.cpp index 7e6900479c..384db8d031 100644 --- a/Numeric/GaussQuadraturePri.cpp +++ b/Numeric/GaussQuadraturePri.cpp @@ -3,54 +3,52 @@ // See the LICENSE.txt file for license information. Please report all // bugs and problems to <gmsh@geuz.org>. +#include "GmshMessage.h" #include "GaussIntegration.h" #include "GaussLegendre1D.h" -#include "stdio.h" IntPt *getGQPriPts(int order); int getNGQPriPts(int order); -IntPt * GQP[] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; +IntPt * GQP[] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, + 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; IntPt *getGQPriPts(int order) -{ +{ int nLin = (order+3)/2; int nTri = getNGQTPts(order); int n = nLin*nTri; int index = order; - if (order >= sizeof(GQP) / sizeof(IntPt*)) + if (order >= (int)(sizeof(GQP) / sizeof(IntPt*))) Msg::Fatal("Increase size of GQP in gauss quadrature prism"); - if(!GQP[index]) - { - double *linPt,*linWt; - IntPt *triPts = getGQTPts(order); -// printf("o = %d n = %d nT = %d nL = %d i = %d\n",order,n,nTri,nLin,index); - gmshGaussLegendre1D(nLin,&linPt,&linWt); - GQP[index] = new IntPt[n]; - int l = 0; - for(int i=0; i < nTri; i++) { - for(int j=0; j < nLin; j++) { - GQP[index][l].pt[0] = triPts[i].pt[0]; - GQP[index][l].pt[1] = triPts[i].pt[1]; - GQP[index][l].pt[2] = linPt[j]; - GQP[index][l++].weight = triPts[i].weight*linWt[j]; -// printf ("%d: %f %f %f %f\n",l-1,triPts[i].pt[0],triPts[i].pt[1],linPt[j],triPts[i].weight*linWt[j]); - } + if(!GQP[index]){ + double *linPt,*linWt; + IntPt *triPts = getGQTPts(order); + gmshGaussLegendre1D(nLin,&linPt,&linWt); + GQP[index] = new IntPt[n]; + int l = 0; + for(int i=0; i < nTri; i++) { + for(int j=0; j < nLin; j++) { + GQP[index][l].pt[0] = triPts[i].pt[0]; + GQP[index][l].pt[1] = triPts[i].pt[1]; + GQP[index][l].pt[2] = linPt[j]; + GQP[index][l++].weight = triPts[i].weight*linWt[j]; } } + } return GQP[index]; } int getNGQPriPts(int order) -{ +{ int nLin = (order+3)/2; int nTri = getNGQTPts(order); - return nLin*nTri; + return nLin * nTri; // if(order == 3)return 8; // if(order == 2)return 8; // if(order < 2) -// return GQPnPt[order]; +// return GQPnPt[order]; // return ((order+3)/2)*((order+3)/2)*((order+3)/2); } diff --git a/Numeric/fullMatrix.h b/Numeric/fullMatrix.h index 0a13ed4cc8..4741dc2c98 100644 --- a/Numeric/fullMatrix.h +++ b/Numeric/fullMatrix.h @@ -541,8 +541,10 @@ class fullMatrix printf("};\n"); } -// specific functions for dgshell - void mult_naiveBlock(const fullMatrix<scalar> &b, const int ncol, const int fcol, const int alpha, const int beta, fullVector<scalar> &c, const int row=0) const + // specific functions for dgshell + void mult_naiveBlock(const fullMatrix<scalar> &b, const int ncol, const int fcol, + const int alpha, const int beta, fullVector<scalar> &c, + const int row=0) const { if(beta != 1) c.scale(beta); @@ -550,7 +552,8 @@ class fullMatrix for(int k = 0; k < _c ; k++) c._data[j] += alpha*(*this)(row, k) * b(k, j); } - void multOnBlock(const fullMatrix<scalar> &b, const int ncol, const int fcol, const int alpha, const int beta, fullVector<scalar> &c) const + void multOnBlock(const fullMatrix<scalar> &b, const int ncol, const int fcol, + const int alpha, const int beta, fullVector<scalar> &c) const #if !defined(HAVE_BLAS) { mult_naiveBlock(b,ncol,fcol,alpha,beta,c); @@ -558,7 +561,8 @@ class fullMatrix #endif ; - void multWithATranspose(const fullVector<scalar> &x, scalar alpha, scalar beta, fullVector<scalar> &y) const + void multWithATranspose(const fullVector<scalar> &x, scalar alpha, scalar beta, + fullVector<scalar> &y) const #if !defined(HAVE_BLAS) { y.scale(beta); @@ -580,7 +584,8 @@ class fullMatrix scalar alpha=1., scalar beta=1.) #if !defined(HAVE_BLAS) { - Msg::Error("gemmWithAtranspose is only available with blas. If blas is not installed please transpose a before used gemm_naive"); + Msg::Error("gemmWithAtranspose is only available with blas. If blas is not " + "installed please transpose a before used gemm_naive"); } #endif ; diff --git a/doc/gmsh.html b/doc/gmsh.html index be84ebc610..cb11d9ee13 100644 --- a/doc/gmsh.html +++ b/doc/gmsh.html @@ -82,7 +82,8 @@ the <a href="/gmsh/doc/LICENSE.txt">GNU General Public License <li>Development version: <ul><li>automated nightly builds (<a href="http://onelab.info/CDash/index.php?project=Gmsh">dashboard</a>): - <a href="/gmsh/bin/Windows/gmsh-svn-Windows.zip">Windows</a>, + Windows <a href="/gmsh/bin/Windows/gmsh-svn-Windows.zip">32 bit</a> / + <a href="/gmsh/bin/Windows/gmsh-svn-Windows64.zip">64 bit</a>, Linux <a href="/gmsh/bin/Linux/gmsh-svn-Linux.tgz">32 bit</a> / <a href="/gmsh/bin/Linux/gmsh-svn-Linux64.tgz">64 bit</a> and -- GitLab