MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
gauss.h
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/libraries/libmbutil/gauss.h,v 1.19 2017/01/12 14:44:04 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 1996-2017
7  *
8  * Pierangelo Masarati <masarati@aero.polimi.it>
9  * Paolo Mantegazza <mantegazza@aero.polimi.it>
10  *
11  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
12  * via La Masa, 34 - 20156 Milano, Italy
13  * http://www.aero.polimi.it
14  *
15  * Changing this copyright notice is forbidden.
16  *
17  * This program is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation (version 2 of the License).
20  *
21  *
22  * This program is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with this program; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  */
31 
32 /* Dati e metodi relativi ai punti di integrazione di Gauss */
33 
34 #ifndef GAUSS_H
35 #define GAUSS_H
36 
37 #include <ac/f2c.h>
38 
39 #include <myassert.h>
40 
41 /* Struttura punto di Gauss.
42  * Ha due campi: la coordinata ed il peso.
43  * Possiede inoltre due metodi per l'accesso in lettura ed un costruttore.
44  * Se si desidera una maggiore robustezza, la si puo' dichiarare come
45  * classe. In tale caso l'accesso in scrittura e' assicurato grazie al
46  * costruttore (e quindi e' riservato a chi genera il punto), mentre la
47  * lettura avviene esclusivamente attraverso i metodi.
48  */
49 
50 struct PntWght {
53 
54  PntWght(doublereal dP, doublereal dW) : dPnt(dP), dWght(dW) { NO_OP; };
55  doublereal dGetPnt(void) const { return dPnt; };
56  doublereal dGetWght(void) const { return dWght; };
57 };
58 
59 
60 /* Versione classe
61 class PntWght {
62  private:
63  doublereal dPnt;
64  doublereal dWght;
65 
66  public:
67  PntWght(doublereal dP, doublereal dW) : dPnt(dP), dWght(dW) { NO_OP; };
68  doublereal dGetPnt(void) const { return dPnt; };
69  doublereal dGetWght(void) const { return dWght; };
70 };
71  */
72 
73 
74 /* Classe set di punti di Gauss.
75  * Viene costruita a partire dal numero di punti che si desidera.
76  * Contiene metodi per la lettura dei singoli punti e pesi, come reali
77  * o come strutture di punti.
78  * Inoltre consente l'accesso ai vettori dei punti (inutile e sconsigliato)
79  *
80  * I dati sono limitati. Sono disponibili per set da 1 a 5 punti.
81  * Se il costruttore viene chiamato con un numero minore di 1,
82  * vengono restituiti i dati relativi ad 1 punto; analogamente,
83  * se il costruttore viene chiamato con un numero maggiore di 5,
84  * vengono restituiti i dati relativi a 5 punti.
85  * Se si utilizza l'iteratore descritto piu' avanti, a parte la perdita
86  * di precisione che ne consegue, l'integrazione non risente quindi
87  * di un valore al di fuori del range previsto.
88  */
89 
90 class NumIntData {
91  protected:
93  public:
94  NumIntData(const integer& i) : iNum(i) {};
95  virtual ~NumIntData(void) {};
96  virtual integer iGetNum(void) const { return iNum; };
97  virtual doublereal dGetPnt(integer i) const = 0;
98  virtual doublereal dGetWght(integer i) const = 0;
99  virtual PntWght Get(integer i) const = 0;
100 };
101 
102 class GaussData : public NumIntData {
103  protected:
104  const doublereal* const pdPnt;
105  const doublereal* const pdWght;
106 
107  public:
108  GaussData(integer iN);
109  doublereal dGetPnt(integer i) const;
110  doublereal dGetWght(integer i) const;
111  PntWght Get(integer i) const;
112  const doublereal* pdGetPnt(void) const;
113  const doublereal* pdGetWght(void) const;
114 };
115 
116 class TrapezoidData : public NumIntData {
117  protected:
118  const doublereal* const pdPnt;
119  const doublereal* const pdWght;
120 
121  public:
123  doublereal dGetPnt(integer i) const;
124  doublereal dGetWght(integer i) const;
125  PntWght Get(integer i) const;
126  const doublereal* pdGetPnt(void) const;
127  const doublereal* pdGetWght(void) const;
128 };
129 
130 
131 /* Classe iteratore sui punti di Gauss.
132  * Dal momento che in ogni caso la classe set di punti di Gauss in ogni caso
133  * e' un'interfaccia per l'accesso diretto ai dati, l'iteratore e' derivato
134  * dalla classe GaussData in contrasto con la filosofia degli iteratori.
135  * In questo modo lo stesso oggetto consente l'accesso diretto ed iterativo
136  * ai dati.
137  * L'utilizzo con la struttura di dati del singolo punto e' consigliato
138  * perche' consente costrutti del tipo:
139  *
140  * GaussDataIterator gdi(n);
141  * PntWght p = gdi.GetFirst();
142  * do {
143  * //... use p.dGetPnt(), p.dGetWght() ...
144  * } while(gdi.fGetNext(p));
145  *
146  * che risultano molto essenziali ed eleganti. E' possibile un uso analogo,
147  * ma meno elegante, dei singoli punti o pesi, qualora fossero richiesti
148  * solo gli uni o gli altri:
149  *
150  * GaussDataIterator gdi(n);
151  * doublereal dp = gdi.dGetFirst();
152  * do {
153  * doublereal dw = gdi.dGetCurrWght();
154  * //...
155  * } while(gdi.fGetNext(dp));
156  *
157  * ATTENZIONE: se si usa dGetCurrPnt() o dGetCurrWght() dopo che fGetNext()
158  * ha ritornato FALSE, il risultato e' unpredictable in quanto il contatore
159  * e' volutamente out of range.
160  */
161 
162 
163 
164 class GaussDataIterator : public GaussData {
165  protected:
167  public:
169  doublereal dGetFirst(integer i = 0) const;
170  PntWght GetFirst(void) const;
171  flag fGetNext(doublereal& d, integer i = 0) const;
172  doublereal dGetCurrPnt(void) const;
173  doublereal dGetCurrWght(void) const;
174  flag fGetNext(PntWght& PW) const;
175 };
176 
177 
179  protected:
182  public:
184  virtual ~NumIntIterator(void);
185  virtual doublereal dGetFirst(integer i = 0) const;
186  virtual PntWght GetFirst(void) const;
187  virtual flag fGetNext(doublereal& d, integer i = 0) const;
188  virtual doublereal dGetCurrPnt(void) const;
189  virtual doublereal dGetCurrWght(void) const;
190  virtual flag fGetNext(PntWght& PW) const;
191 };
192 
193 #endif /* GAUSS_H */
const doublereal *const pdWght
Definition: gauss.h:119
long int flag
Definition: mbdyn.h:43
virtual doublereal dGetPnt(integer i) const =0
flag fGetNext(doublereal &d, integer i=0) const
Definition: gauss.cc:205
virtual PntWght GetFirst(void) const
Definition: gauss.cc:283
PntWght GetFirst(void) const
Definition: gauss.cc:198
NumIntData & data
Definition: gauss.h:181
doublereal dGetWght(integer i) const
Definition: gauss.cc:103
const doublereal * pdGetPnt(void) const
Definition: gauss.cc:162
virtual PntWght Get(integer i) const =0
NumIntData(const integer &i)
Definition: gauss.h:94
#define NO_OP
Definition: myassert.h:74
PntWght Get(integer i) const
Definition: gauss.cc:155
virtual flag fGetNext(doublereal &d, integer i=0) const
Definition: gauss.cc:290
const doublereal * pdGetPnt(void) const
Definition: gauss.cc:117
virtual ~NumIntIterator(void)
Definition: gauss.cc:264
doublereal dGetFirst(integer i=0) const
Definition: gauss.cc:185
doublereal dGetWght(void) const
Definition: gauss.h:56
const doublereal *const pdPnt
Definition: gauss.h:118
const doublereal * pdGetWght(void) const
Definition: gauss.cc:123
doublereal dGetWght(integer i) const
Definition: gauss.cc:148
GaussData(integer iN)
Definition: gauss.cc:87
doublereal dPnt
Definition: gauss.h:51
const doublereal * pdGetWght(void) const
Definition: gauss.cc:168
PntWght Get(integer i) const
Definition: gauss.cc:110
doublereal dGetCurrWght(void) const
Definition: gauss.cc:232
integer iCurr
Definition: gauss.h:180
GaussDataIterator(integer iN)
Definition: gauss.cc:178
doublereal dGetPnt(integer i) const
Definition: gauss.cc:141
PntWght(doublereal dP, doublereal dW)
Definition: gauss.h:54
doublereal dGetCurrPnt(void) const
Definition: gauss.cc:224
TrapezoidData(integer iN)
Definition: gauss.cc:133
doublereal dGetPnt(void) const
Definition: gauss.h:55
integer iCurr
Definition: gauss.h:166
virtual integer iGetNum(void) const
Definition: gauss.h:96
NumIntIterator(NumIntData &d)
Definition: gauss.cc:257
doublereal dGetPnt(integer i) const
Definition: gauss.cc:96
virtual doublereal dGetFirst(integer i=0) const
Definition: gauss.cc:270
virtual doublereal dGetCurrPnt(void) const
Definition: gauss.cc:309
double doublereal
Definition: colamd.c:52
doublereal dWght
Definition: gauss.h:52
long int integer
Definition: colamd.c:51
const doublereal *const pdWght
Definition: gauss.h:105
const doublereal *const pdPnt
Definition: gauss.h:104
virtual doublereal dGetCurrWght(void) const
Definition: gauss.cc:317
integer iNum
Definition: gauss.h:92
virtual doublereal dGetWght(integer i) const =0
Definition: gauss.h:50
virtual ~NumIntData(void)
Definition: gauss.h:95