\psfrag{T4}{\scriptsize $T_4$}
\begin{document}
-% \todo{
-\input{titelseite}
-% }
+\todo{
+% \input{titelseite}
+}
\tableofcontents
\clearpage
In dieser Arbeit beschäftigen wir uns mit der Randelementmethode für die homogene Laplace-Gleichung mit Dirichlet-Randbedingungen
\begin{align}
- \varDelta u &= 0 \text{ in } \Omega \subset \R^3,\\
- u &= f \text{ auf } \Gamma,
+ u &= g \text{ auf } \Gamma,
\end{align}
wobei $\varDelta u := \partial_x^2u+\partial_y^2u$ den Laplace-Operator bezeichnet und $\Omega \subset \R^3$ eine beschränkte Teilmenge von $\R^3$ mit Lipschitz-Rand $\Gamma := \partial \Omega$ ist.\\
\todo{
\begin{itemize}
\item kuerzer als in Ferraz-DA \cite{fer:errbem}
- \item insbesondere fastregulaere Triangulierung + K-Gitter (vgl. Ferraz-DA)
+ \item insbesondere fastregulaere Triangulierung + K-Gitter (vgl. Ferraz-DA \cite{fer:errbem})
\end{itemize}
}
\subsection{Netze}
Für die Diskretisierung des Problems wollen wir nun einige Begriffe definieren.
\begin{defi}[Rechteck]
- Wir nennen $T \subseteq \R^3$ ein achsenorientiertes Rechteck, wenn sich aus einem Eckpunkt $C \in \R^3$, alle weiteren Knoten schreiben lassen als $\{ C + a\cdot A, C + a \cdot A + b \cdot B , C + b \cdot B\}$, wobei Vektoren $A,B \in \{(1,0,0)^T,(0,1,0)^T,(0,0,1)^T\}$. und Skalare $a,b\neq0$ sind und der Flächeninhalt $|T| > 0$. \todo{Kanten/Seiten definieren? Anmerken, dass ein Kante als Vektor nur eine dim $\neq 0$?}\\Weiterhin sei
- \begin{align}
- C_T &:= \{ c_i:1\leq i \leq 4\}
- \end{align}
-die Menge der Eckpunkte des Rechtecks.
+ Wir nennen $T \subseteq \R^3$ ein achsenorientiertes Rechteck, wenn sich aus einem Eckpunkt $C \in \R^3$, alle weiteren Knoten schreiben lassen als $\{ C + a\cdot A, C + a \cdot A + b \cdot B , C + b \cdot B\}$, wobei Vektoren $A,B \in \{(1,0,0)^T,(0,1,0)^T,(0,0,1)^T\}$. und Skalare $a,b\neq0$ und der Flächeninhalt $|T| > 0$. \todo{Kanten/Seiten definieren? Anmerken, dass ein Kante als Vektor nur eine dim $\neq 0$?}
\end{defi}
\begin{defi}[Partition]
\begin{itemize}
\item $\overline{\Gamma} = \bigcup_{j=1}^NT_j$
\item alle Elemente aus $T_{\ell}$ abgeschlossene achsenorientierte Rechtecke
- \item $T_j \cap T_k$ ist entweder ein Knoten, eine Seite oder keines von beiden mit $T_j,T_k\in\T_{\ell}$ und $T_j\neq T_k$
- \item auf einer Kante liegt maximal ein hängender Knoten
+ \item $\abs{T_j \cap T_k}=0$ für den Schnitt zweier Elemente $T_j,T_k\in\T_{\ell}$ mit $T_j\neq T_k$
\end{itemize}
+Hier bezeichne $\abs{\cdot}$ das 2-Dimensionale Oberflächenmaß für Rechtecke.
\end{defi}
-Daraus folgt, dass pro Seite eines Rechtecks
+Weiterhin wollen wir in dieser Arbeit Aufgrund der Netzstabilität noch fordern, dass der Schnitt $T_j \cap T_k$ zweier Elemente $T_j,T_k\in\T_{\ell}$ mit $T_j\neq T_k$
\begin{itemize}
- \item entweder kein Nachbarelement
- \item oder genau ein Nachbarelement
- \item oder genau zwei Nachbarelemente
+ \item leer ist,
+ \item oder Knoten von $T_j$ und $T_k$,
+ \item oder Kante von $T_j$ und $T_k$,
+ \item oder o.B.d.A. Kante von $T_j$ und Kantenhälfte von $T_k$.
\end{itemize}
-anliegen können. Weiterhin wollen wir auch verlangen, dass der auf der Kante liegende Knoten die Seite halbiert, die beiden Nachbarlemente also die gleiche Seitenlänge an der anliegenden Kante haben.
+
\begin{figure}[ht]
\centering
\end{figure}
\begin{defi}
- Für eine Partition $\T_{\ell}$ bezeichnen wir mit $\hat \T_{\ell}$ jene Partition die entsteht, wenn alle Elemente isotrop verfeinert werden.
+ Für eine Partition $\T_{\ell}$ bezeichnen wir mit $\widehat \T_{\ell}$ jene Partition die entsteht, wenn alle Elemente isotrop verfeinert werden.
\end{defi}
\begin{alg}[Verfeinern] \label{alg:refine} Sei $\T_{\ell} = \{T_1,T_2,\dots,T_N\}$ eine Partition und $marked$ eine gegebene Markierung wobei $marked_j \hat =$ Markierung von $T_j$. Gehe so vor:
\begin{enumerate}
\subsection{Gauss-Quadratur}
\begin{displaymath}
- \int_0^l f(x) dx \approx \sum_{k=0}^n w_kf(x_k) =: \Q(f)
+ \int_0^s f(x) dx \approx \sum_{k=0}^n w_kf(x_k) =: \Q(f)
\end{displaymath}
\subsection{Quadratur über ein Element}
}
\subsection{Fehlerschätzer}
-In diesem Abschnitt definieren wir die a-posteriori Fehlerschätzer, die wir im Folgenden zur Steuerung des adaptiven Algorithmus einsetzen werden. Wir verwenden dazu die $h-h/2$ Strategie aus Ferraz-Leite (Provet. 108),wo die folgende Aussage bewiesen wird.
+In diesem Abschnitt definieren wir die a-posteriori Fehlerschätzer, die wir im Folgenden zur Steuerung des adaptiven Algorithmus einsetzen werden. Wir verwenden dazu die $h-h/2$ Strategie aus Ferraz-Leite, wo die folgende Aussage bewiesen wird.
\begin{defi}Es bezeichne $\phi$ die Lösung von Formel , $\phi_{\ell}$ die Galerkinlösung auf dem Gitter $\T_{\ell}$ und $\hat \phi_{\ell}$ die Galerkinlösung auf dem uniform verfeinerten Gitter $\hat \T_{\ell}$. Dann gilt, der Schätzer
\begin{align}
\renewcommand{\ttdefault}{pcr}
\definecolor{lst_green}{rgb}{0,0.6,0}
-\definecolor{lst_lgreen}{rgb}{0.6,0.6,0}
+\definecolor{lst_dbleen}{rgb}{0,0.3,0.3}
\definecolor{lst_gray}{rgb}{0.6,0.6,0.6}
+\definecolor{lst_dgray}{rgb}{0.6,0.6,0.6}
\definecolor{lst_mauve}{rgb}{0.6,0,0.6}
\definecolor{lst_dmauve}{rgb}{0.3,0,0.3}
\definecolor{lst_blue}{rgb}{0,0,0.6}
% breakatwhitespace=false,
% title=\lstname,
% identifierstyle=\color{red}, %
+% morekeywords=[4]{int,double,void}, %Typen
keywordstyle=[1]\color{lst_mauve}\bfseries, %buildin commands
keywordstyle=[2]\bfseries\color{lst_dmauve}, %lib commands
keywordstyle=[3]\color{lst_bleen}, %globals
- keywordstyle=[4]\bfseries, %Functions
- commentstyle=\color{lst_green},
+ keywordstyle=[4]\bfseries\color{lst_dbleen}, %Typen
+ keywordstyle=[5]\bfseries, %Functions
+ commentstyle=\color{lst_dgray},
stringstyle=\color{lst_blue},
% escapeinside={\%*}{*)},
-extendedchars=true
+ extendedchars=true
}
\vspace{1.5cm}
unter der Anleitung von \\
Ao.Univ.Prof. Dipl.Math. Dr.techn. Dirk Praetorius \\
-Dipl.-Ing. Michael Karkulik\\
+Projektass.(FWF) Dipl.-Ing. Michael Karkulik\\
\vspace{1cm}
durch \\
Peter Matthias Schaefer \\
Schopenhauerstr. 48/2 \\
1180 Wien
-
\end{center}
\end{titlepage}
function str = export_gauss(start,stop,steps,file)
str = ['#ifndef GAUSS_NODES' char(10) ...
-'#define GAUSS_NODES' char(10) ...
+... '#define GAUSS_NODES' char(10) ...
'/*' char(10) ...
-' * Gaus Auswertungspunkte auf dem Intervall [0 1]' char(10) ...
+' * Gaus Auswertungspunkte auf dem Intervall'...
+ ' [' num2str(start) ' ' num2str(stop) ']' char(10) ...
' *' char(10) ...
' * 2^$NUM entspricht der Anzahl der Auswertungspunkte' char(10) ...
' *' char(10) ...
-' * Anzahl der Punkte: gauss_size[$NUM]' char(10) ...
-' * Auswertungstellen: gauss_nodes[$NUM][$I].n ($I aktuelle Position)' char(10) ...
-' * Gewichte: gauss_nodes[$NUM][$I].w ($I aktuelle Position)' char(10) ...
+' * Anzahl der Punkte: GAUSS_SIZE[$NUM]' char(10) ...
+' * Auswertungstellen: GAUSS_NODES[$NUM][$I].n ($I aktuelle Position)' char(10) ...
+' * Gewichte: GAUSS_NODES[$NUM][$I].w ($I aktuelle Position)' char(10) ...
' *' char(10) ...
' * Beispiel:' char(10) ...
' *' char(10) ...
' * > double sol;' char(10) ...
-' * > for(int i=0;i<gauss_size[3];++i)' char(10) ...
-' * > sol += sin(gauss_nodes[3][i].n) * gauss_nodes[3][i].w;' char(10) ...
+' * > for(int i=0;i<GAUSS_SIZE[3];++i)' char(10) ...
+' * > sol += sin(GAUSS_NODES[3][i].n) * GAUSS_NODES[3][i].w;' char(10) ...
' *' char(10) ...
-' * Peter Schaefer (Automatisch generiert durch export_gaus.m )' char(10) ...
+' * Peter Schaefer (Automatisch generiert durch export_gauss.m )' char(10) ...
' */' char(10) char(10) ...
-'typedef struct _gauss{' char(10) ...
+'typedef struct _gauss {' char(10) ...
' double n;' char(10) ...
' double w;' char(10) ...
'} gauss;' char(10) ...
preci = 21;
-str2 =['gauss gauss_nodes[][' num2str(max(steps)) '] = {' char(10)];
-str1 =['double gauss_size [] = {'];
+str2 =['gauss const GAUSS_NODES[][' num2str(max(steps)) '] = {' char(10)];
+str1 =['double const GAUSS_SIZE[] = { '];
for i = steps(1:end)
- str1 = [str1 num2str(i) ','];
+ str1 = [str1 num2str(i) ', '];
str2 = [str2 '{' char(10)];
[n w] = gauss(i,start,stop);
for j = 1:length(n)
- str2 = [str2 char(9) '{' num2str(n(j),preci) ',' num2str(w(j),preci) '},' char(10)];
+ str2 = [str2 char(9) '{ ' num2str(n(j),preci) ', ' num2str(w(j),preci) ' },' char(10)];
end
str2 = [str2(1:end-2) '},' char(10)];
end
-str1 = [str1(1:end-1) '};' char(10)];
+str1 = [str1(1:end-2) ' };' char(10)];
str2 = [str2(1:end-2) '};' char(10)];
str = [str str1 char(10) str2 char(10) '#endif' char(10)];
#ifndef GAUSS_NODES
-#define GAUSS_NODES
/*
* Gaus Auswertungspunkte auf dem Intervall [0 1]
*
* 2^$NUM entspricht der Anzahl der Auswertungspunkte
*
- * Anzahl der Punkte: gauss_size[$NUM]
- * Auswertungstellen: gauss_nodes[$NUM][$I].n ($I aktuelle Position)
- * Gewichte: gauss_nodes[$NUM][$I].w ($I aktuelle Position)
+ * Anzahl der Punkte: GAUSS_SIZE[$NUM]
+ * Auswertungstellen: GAUSS_NODES[$NUM][$I].n ($I aktuelle Position)
+ * Gewichte: GAUSS_NODES[$NUM][$I].w ($I aktuelle Position)
*
* Beispiel:
*
* > double sol;
- * > for(int i=0;i<gauss_size[3];++i)
- * > sol += sin(gauss_nodes[3][i].n) * gauss_nodes[3][i].w;
+ * > for(int i=0;i<GAUSS_SIZE[3];++i)
+ * > sol += sin(GAUSS_NODES[3][i].n) * GAUSS_NODES[3][i].w;
*
- * Peter Schaefer (Automatisch generiert durch export_gaus.m )
+ * Peter Schaefer (Automatisch generiert durch export_gauss.m )
*/
-typedef struct _gauss{
+typedef struct _gauss {
double n;
double w;
} gauss;
-double gauss_size [] = {2,4,8,16,32};
+double const GAUSS_SIZE[] = { 1, 2, 4, 8, 16, 32 };
-gauss gauss_nodes[][32] = {
+gauss const GAUSS_NODES[][32] = {
{
- {0.211324865405187078959,0.499999999999999888978},
- {0.788675134594812865529,0.499999999999999888978}},
+ { 0.5, 1 }},
{
- {0.0694318442029737137311,0.173927422568727063634},
- {0.330009478207571815833,0.32607257743127315841},
- {0.669990521792428128656,0.326072577431272991877},
- {0.93056815579702634178,0.173927422568727008123}},
+ { 0.211324865405187078959, 0.499999999999999888978 },
+ { 0.788675134594812865529, 0.499999999999999888978 }},
{
- {0.0198550717512319119251,0.0506142681451882195387},
- {0.101666761293186636017,0.111190517226687282659},
- {0.237233795041835615613,0.15685332293894366229},
- {0.408282678752175054449,0.181341891689181078373},
- {0.591717321247825056574,0.18134189168918116164},
- {0.762766204958164384387,0.156853322938943606779},
- {0.898333238706813475005,0.111190517226687282659},
- {0.980144928248767977053,0.0506142681451882542332}},
+ { 0.0694318442029737137311, 0.173927422568727063634 },
+ { 0.330009478207571815833, 0.32607257743127315841 },
+ { 0.669990521792428128656, 0.326072577431272991877 },
+ { 0.93056815579702634178, 0.173927422568727008123 }},
{
- {0.00529953250417525278948,0.0135762297058772823943},
- {0.0277124884633836443548,0.0311267619693235096656},
- {0.0671843988060840668908,0.0475792558412463095774},
- {0.12229779582249833414,0.0623144856277671119194},
- {0.191061877798677837159,0.0747979944082881875733},
- {0.270991611171386259649,0.0845782596975013928331},
- {0.359198224610370597798,0.0913017075224623053664},
- {0.452493745081181397705,0.0947253052275341400623},
- {0.547506254918818657806,0.0947253052275344453736},
- {0.640801775389629457713,0.0913017075224619306661},
- {0.729008388828613518307,0.0845782596975014205887},
- {0.808938122201321996307,0.0747979944082887010515},
- {0.877702204177501554838,0.0623144856277674172307},
- {0.932815601193915933109,0.0475792558412462263107},
- {0.972287511536616300134,0.0311267619693240300827},
- {0.994700467495824858233,0.0135762297058771956582}},
+ { 0.0198550717512319119251, 0.0506142681451882195387 },
+ { 0.101666761293186636017, 0.111190517226687282659 },
+ { 0.237233795041835615613, 0.15685332293894366229 },
+ { 0.408282678752175054449, 0.181341891689181078373 },
+ { 0.591717321247825056574, 0.18134189168918116164 },
+ { 0.762766204958164384387, 0.156853322938943606779 },
+ { 0.898333238706813475005, 0.111190517226687282659 },
+ { 0.980144928248767977053, 0.0506142681451882542332 }},
{
- {0.00136806907525910403933,0.00350930500473500784839},
- {0.00719424422736580915227,0.00813719736545263742922},
- {0.0176188722062468050567,0.0126960326546312531754},
- {0.0325469620311301111037,0.0171369314565104867432},
- {0.0518394221169737878796,0.0214179490111136711095},
- {0.0753161931337150702959,0.0254990296311881116387},
- {0.102758102016028862735,0.0293420467392677027096},
- {0.133908940629855199855,0.0329111113881807790249},
- {0.168477866534892439798,0.0361728970544242522944},
- {0.206142121379618625809,0.0390969478935349543103},
- {0.246550045533885375804,0.041655962113473471442},
- {0.289324361934682361408,0.0438260465022020373471},
- {0.334065698858936110938,0.0455869393478821535726},
- {0.380356318873931620317,0.0469221995404019084908},
- {0.427764019208601742328,0.0478193600396371459871},
- {0.475846167156130650522,0.0482700442573640656208},
- {0.524153832843869071922,0.0482700442573634549981},
- {0.572235980791398368694,0.0478193600396370696592},
- {0.619643681126068490705,0.046922199540402220741},
- {0.66593430114106377804,0.0455869393478817580556},
- {0.710675638065317860637,0.0438260465022019540804},
- {0.753449954466114735219,0.041655962113473277153},
- {0.793857878620381374191,0.0390969478935357661609},
- {0.831522133465107504691,0.0361728970544242176},
- {0.866091059370144966678,0.032911111388180640247},
- {0.897241897983971248287,0.029342046739267681893},
- {0.924683806866285151749,0.0254990296311881359248},
- {0.94816057788302621212,0.021417949011113403962},
- {0.967453037968869833385,0.0171369314565111112436},
- {0.982381127793753194943,0.0126960326546313780754},
- {0.992805755772633968803,0.00813719736545309192677},
- {0.998631930924740673916,0.00350930500473496925079}}};
+ { 0.00529953250417525278948, 0.0135762297058772823943 },
+ { 0.0277124884633836443548, 0.0311267619693235096656 },
+ { 0.0671843988060840668908, 0.0475792558412463095774 },
+ { 0.12229779582249833414, 0.0623144856277671119194 },
+ { 0.191061877798677837159, 0.0747979944082881875733 },
+ { 0.270991611171386259649, 0.0845782596975013928331 },
+ { 0.359198224610370597798, 0.0913017075224623053664 },
+ { 0.452493745081181397705, 0.0947253052275341400623 },
+ { 0.547506254918818657806, 0.0947253052275344453736 },
+ { 0.640801775389629457713, 0.0913017075224619306661 },
+ { 0.729008388828613518307, 0.0845782596975014205887 },
+ { 0.808938122201321996307, 0.0747979944082887010515 },
+ { 0.877702204177501554838, 0.0623144856277674172307 },
+ { 0.932815601193915933109, 0.0475792558412462263107 },
+ { 0.972287511536616300134, 0.0311267619693240300827 },
+ { 0.994700467495824858233, 0.0135762297058771956582 }},
+{
+ { 0.00136806907525910403933, 0.00350930500473500784839 },
+ { 0.00719424422736580915227, 0.00813719736545263742922 },
+ { 0.0176188722062468050567, 0.0126960326546312531754 },
+ { 0.0325469620311301111037, 0.0171369314565104867432 },
+ { 0.0518394221169737878796, 0.0214179490111136711095 },
+ { 0.0753161931337150702959, 0.0254990296311881116387 },
+ { 0.102758102016028862735, 0.0293420467392677027096 },
+ { 0.133908940629855199855, 0.0329111113881807790249 },
+ { 0.168477866534892439798, 0.0361728970544242522944 },
+ { 0.206142121379618625809, 0.0390969478935349543103 },
+ { 0.246550045533885375804, 0.041655962113473471442 },
+ { 0.289324361934682361408, 0.0438260465022020373471 },
+ { 0.334065698858936110938, 0.0455869393478821535726 },
+ { 0.380356318873931620317, 0.0469221995404019084908 },
+ { 0.427764019208601742328, 0.0478193600396371459871 },
+ { 0.475846167156130650522, 0.0482700442573640656208 },
+ { 0.524153832843869071922, 0.0482700442573634549981 },
+ { 0.572235980791398368694, 0.0478193600396370696592 },
+ { 0.619643681126068490705, 0.046922199540402220741 },
+ { 0.66593430114106377804, 0.0455869393478817580556 },
+ { 0.710675638065317860637, 0.0438260465022019540804 },
+ { 0.753449954466114735219, 0.041655962113473277153 },
+ { 0.793857878620381374191, 0.0390969478935357661609 },
+ { 0.831522133465107504691, 0.0361728970544242176 },
+ { 0.866091059370144966678, 0.032911111388180640247 },
+ { 0.897241897983971248287, 0.029342046739267681893 },
+ { 0.924683806866285151749, 0.0254990296311881359248 },
+ { 0.94816057788302621212, 0.021417949011113403962 },
+ { 0.967453037968869833385, 0.0171369314565111112436 },
+ { 0.982381127793753194943, 0.0126960326546313780754 },
+ { 0.992805755772633968803, 0.00813719736545309192677 },
+ { 0.998631930924740673916, 0.00350930500473496925079 }}};
#endif
double * E = mxGetPr(prhs[1]);\r
\r
int type = (int) *(mxGetPr(prhs[3]));\r
- double * zeta ={};\r
+ double * zeta = {0};\r
\r
#if DEBUG >= 4\r
mexPrintf("ZetaM:%d\nZetaN:%d\n",mxGetM(prhs[2]),mxGetN(prhs[2]));\r
#endif\r
- if(mxGetN(prhs[2])>1)\r
- zeta = mxGetPr(prhs[2])+1;\r
+ if (mxGetN(prhs[2]) > 1)\r
+ zeta = mxGetPr(prhs[2]) + 1;\r
\r
int quad = (int) *mxGetPr(prhs[2]);\r
quad = setQuad(quad);\r
mexPrintf("Berechnungstyp:%d\nQuadpunkte:%d\n",type,quad);\r
#endif\r
\r
-\r
// Welche Funktion soll verwendet werden\r
double (*ctypeP)(double, double, double, double, double, double, double,\r
double*);\r
#include "gauss.hpp"\r
\r
//#define quad 2 // Anzahl der Quadratur Auswertungstellen 2^quad\r
-int quad=2;\r
+int quad = 2;\r
//using namespace std;\r
\r
-int setQuad(int q){\r
+int setQuad(int q) {\r
quad = q;\r
- return gauss_size[quad];\r
+ return GAUSS_SIZE[quad];\r
}\r
\r
int sign(double);\r
\r
double sol = 0;\r
\r
- for (int i = 0; i < gauss_size[quad]; ++i) {\r
- sol += pow(\r
- (x - gauss_nodes[quad][i].n * y)\r
- * (x - gauss_nodes[quad][i].n * y) + l * l, p)\r
- * gauss_nodes[quad][i].w * y;\r
- }\r
+ for (int i = 0; i < GAUSS_SIZE[quad]; ++i) {\r
+sol += pow(\r
+ (x - GAUSS_NODES[quad][i].n * y)\r
+ * (x - GAUSS_NODES[quad][i].n * y) + l * l, p)\r
+ * GAUSS_NODES[quad][i].w * y;\r
+}\r
\r
return sol;\r
}\r
double l) {\r
\r
double sol = 0;\r
- for (int i = 0; i < gauss_size[quad]; ++i) {\r
- for (int j = 0; j < gauss_size[quad]; ++j) {\r
- sol += pow(\r
- (x1 - y1 * gauss_nodes[quad][i].n)\r
- * (x1 - y1 * gauss_nodes[quad][i].n)\r
- + (x2 - y2 * gauss_nodes[quad][j].n)\r
- * (x2 - y2 * gauss_nodes[quad][j].n)\r
- + l * l, p) * y1 * gauss_nodes[quad][i].w * y2\r
- * gauss_nodes[quad][j].w;\r
- }\r
+ for (int i = 0; i < GAUSS_SIZE[quad]; ++i) {\r
+ for (int j = 0; j < GAUSS_SIZE[quad]; ++j) {\r
+sol += pow(\r
+ (x1 - y1 * GAUSS_NODES[quad][i].n)\r
+ * (x1 - y1 * GAUSS_NODES[quad][i].n)\r
+ + (x2 - y2 * GAUSS_NODES[quad][j].n)\r
+ * (x2 - y2 * GAUSS_NODES[quad][j].n)\r
+ + l * l, p) * y1 * GAUSS_NODES[quad][i].w * y2\r
+ * GAUSS_NODES[quad][j].w;\r
}\r
+}\r
return sol;\r
}\r
\r
//Fall 2\r
double sol = 0;\r
\r
- for (int i = 0; i < gauss_size[quad]; ++i) {\r
- for (int j = 0; j < gauss_size[quad]; ++j) {\r
- sol += G_AY1Y2(-0.5, t + d1, v + d2, b * gauss_nodes[quad][i].n,\r
- d * gauss_nodes[quad][j].n, d3) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol -= G_AY1Y2(-0.5, d1, v + d2, b * gauss_nodes[quad][i].n,\r
- d * gauss_nodes[quad][j].n, d3) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol -= G_AY1Y2(-0.5, t + d1, d2, b * gauss_nodes[quad][i].n,\r
- d * gauss_nodes[quad][j].n, d3) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol += G_AY1Y2(-0.5, d1, d2, b * gauss_nodes[quad][i].n,\r
- d * gauss_nodes[quad][j].n, d3) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- }\r
+ for (int i = 0; i < GAUSS_SIZE[quad]; ++i) {\r
+ for (int j = 0; j < GAUSS_SIZE[quad]; ++j) {\r
+sol += G_AY1Y2(-0.5, t + d1, v + d2, b * GAUSS_NODES[quad][i].n,\r
+ d * GAUSS_NODES[quad][j].n, d3) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol -= G_AY1Y2(-0.5, d1, v + d2, b * GAUSS_NODES[quad][i].n,\r
+ d * GAUSS_NODES[quad][j].n, d3) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol -= G_AY1Y2(-0.5, t + d1, d2, b * GAUSS_NODES[quad][i].n,\r
+ d * GAUSS_NODES[quad][j].n, d3) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol += G_AY1Y2(-0.5, d1, d2, b * GAUSS_NODES[quad][i].n,\r
+ d * GAUSS_NODES[quad][j].n, d3) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
}\r
+}\r
\r
return sol;\r
\r
//Fall 3\r
double sol = 0;\r
\r
- for (int i = 0; i < gauss_size[quad]; ++i) {\r
- for (int j = 0; j < gauss_size[quad]; ++j) {\r
-\r
- sol += G_AY2X2(t * gauss_nodes[quad][j].n + d1, v + d2,\r
- b * gauss_nodes[quad][i].n, d, d3) * t * b\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol -= G_AY2X2(t * gauss_nodes[quad][j].n + d1, v + d2,\r
- b * gauss_nodes[quad][i].n, 0, d3) * t * b\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol -= G_AY2X2(t * gauss_nodes[quad][j].n + d1, d2,\r
- b * gauss_nodes[quad][i].n, d, d3) * t * b\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol += G_AY2X2(t * gauss_nodes[quad][j].n + d1, d2,\r
- b * gauss_nodes[quad][i].n, 0, d3) * t * b\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
+ for (int i = 0; i < GAUSS_SIZE[quad]; ++i) {\r
+ for (int j = 0; j < GAUSS_SIZE[quad]; ++j) {\r
+\r
+ sol += G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, v + d2,\r
+ b * GAUSS_NODES[quad][i].n, d, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol -= G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, v + d2,\r
+ b * GAUSS_NODES[quad][i].n, 0, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol -= G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, d2,\r
+ b * GAUSS_NODES[quad][i].n, d, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol += G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, d2,\r
+ b * GAUSS_NODES[quad][i].n, 0, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
\r
if (sol\r
!= sol /*|| fabs(sol) == numeric_limits<double>::infinity()*/) {\r
/* cout << "->(" << i << "," << j << ")" << endl;\r
cout << "-|("\r
- << G_AY2X2(t * gauss_nodes[quad][j].n + d1, v + d2,\r
- b * gauss_nodes[quad][i].n, d, d3) * t * b\r
- * gauss_nodes[quad][i].w\r
- * gauss_nodes[quad][j].w << ","\r
- << G_AY2X2(t * gauss_nodes[quad][j].n + d1, v + d2,\r
- b * gauss_nodes[quad][i].n, 0, d3) * t * b\r
- * gauss_nodes[quad][i].w\r
- * gauss_nodes[quad][j].w << ","\r
- << G_AY2X2(t * gauss_nodes[quad][j].n + d1, d2,\r
- b * gauss_nodes[quad][i].n, d, d3) * t * b\r
- * gauss_nodes[quad][i].w\r
- * gauss_nodes[quad][j].w << ","\r
- << G_AY2X2(t * gauss_nodes[quad][j].n + d1, d2,\r
- b * gauss_nodes[quad][i].n, 0, d3) * t * b\r
- * gauss_nodes[quad][i].w\r
- * gauss_nodes[quad][j].w << ")" << endl;\r
+ << G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, v + d2,\r
+ b * GAUSS_NODES[quad][i].n, d, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w\r
+ * GAUSS_NODES[quad][j].w << ","\r
+ << G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, v + d2,\r
+ b * GAUSS_NODES[quad][i].n, 0, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w\r
+ * GAUSS_NODES[quad][j].w << ","\r
+ << G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, d2,\r
+ b * GAUSS_NODES[quad][i].n, d, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w\r
+ * GAUSS_NODES[quad][j].w << ","\r
+ << G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, d2,\r
+ b * GAUSS_NODES[quad][i].n, 0, d3) * t * b\r
+ * GAUSS_NODES[quad][i].w\r
+ * GAUSS_NODES[quad][j].w << ")" << endl;\r
cout << "||"\r
- << t * b * gauss_nodes[quad][i].w\r
- * gauss_nodes[quad][j].w\r
+ << t * b * GAUSS_NODES[quad][i].w\r
+ * GAUSS_NODES[quad][j].w\r
<< "||-----------------" << endl;\r
\r
cout << "? ("\r
- << G_AY2X2(t * gauss_nodes[quad][j].n + d1, d2,\r
- b * gauss_nodes[quad][i].n, 0, d3) << ","\r
- << t * b * gauss_nodes[quad][i].w\r
- * gauss_nodes[quad][j].w << ")" << endl;\r
+ << G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, d2,\r
+ b * GAUSS_NODES[quad][i].n, 0, d3) << ","\r
+ << t * b * GAUSS_NODES[quad][i].w\r
+ * GAUSS_NODES[quad][j].w << ")" << endl;\r
\r
cout << endl;*/\r
return sol;\r
double sol = 0;\r
int i, j, k, l;\r
\r
- for (i = 0; i < gauss_size[quad]; ++i) {\r
- for (j = 0; j < gauss_size[quad]; ++j) {\r
- for (k = 0; k < gauss_size[quad]; ++k) {\r
- for (l = 0; l < gauss_size[quad]; ++l) {\r
- sol += f_A(b * gauss_nodes[quad][i].n,\r
- d * gauss_nodes[quad][j].n,\r
- (t) * gauss_nodes[quad][k].n + d1,\r
- (v) * gauss_nodes[quad][l].n + d2, d3) * b\r
- * gauss_nodes[quad][i].w * d\r
- * gauss_nodes[quad][j].w * t\r
- * gauss_nodes[quad][k].w * v\r
- * gauss_nodes[quad][l].w;\r
- }\r
+ for (i = 0; i < GAUSS_SIZE[quad]; ++i) {\r
+ for (j = 0; j < GAUSS_SIZE[quad]; ++j) {\r
+ for (k = 0; k < GAUSS_SIZE[quad]; ++k) {\r
+ for (l = 0; l < GAUSS_SIZE[quad]; ++l) {\r
+sol += f_A(b * GAUSS_NODES[quad][i].n,\r
+ d * GAUSS_NODES[quad][j].n,\r
+ (t) * GAUSS_NODES[quad][k].n + d1,\r
+ (v) * GAUSS_NODES[quad][l].n + d2, d3) * b\r
+ * GAUSS_NODES[quad][i].w * d\r
+ * GAUSS_NODES[quad][j].w * t\r
+ * GAUSS_NODES[quad][k].w * v\r
+ * GAUSS_NODES[quad][l].w;\r
}\r
}\r
}\r
+}\r
\r
return sol;\r
}\r
//Fall 2\r
double sol = 0;\r
\r
- for (int i = 0; i < gauss_size[quad]; ++i) {\r
- for (int j = 0; j < gauss_size[quad]; ++j) {\r
- sol += G_AY1Y2(-0.5, t + d2, v + d3, d * gauss_nodes[quad][j].n, 0,\r
- b * gauss_nodes[quad][i].n - d1) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol -= G_AY1Y2(-0.5, d2, v + d3, d * gauss_nodes[quad][j].n, 0,\r
- b * gauss_nodes[quad][i].n - d1) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol -= G_AY1Y2(-0.5, t + d2, d3, d * gauss_nodes[quad][j].n, 0,\r
- b * gauss_nodes[quad][i].n - d1) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- sol += G_AY1Y2(-0.5, d2, d3, d * gauss_nodes[quad][j].n, 0,\r
- b * gauss_nodes[quad][i].n - d1) * b * d\r
- * gauss_nodes[quad][i].w * gauss_nodes[quad][j].w;\r
- }\r
+ for (int i = 0; i < GAUSS_SIZE[quad]; ++i) {\r
+ for (int j = 0; j < GAUSS_SIZE[quad]; ++j) {\r
+sol += G_AY1Y2(-0.5, t + d2, v + d3, d * GAUSS_NODES[quad][j].n, 0,\r
+ b * GAUSS_NODES[quad][i].n - d1) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol -= G_AY1Y2(-0.5, d2, v + d3, d * GAUSS_NODES[quad][j].n, 0,\r
+ b * GAUSS_NODES[quad][i].n - d1) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol -= G_AY1Y2(-0.5, t + d2, d3, d * GAUSS_NODES[quad][j].n, 0,\r
+ b * GAUSS_NODES[quad][i].n - d1) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
+ sol += G_AY1Y2(-0.5, d2, d3, d * GAUSS_NODES[quad][j].n, 0,\r
+ b * GAUSS_NODES[quad][i].n - d1) * b * d\r
+ * GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;\r
}\r
+}\r
\r
return sol;\r
\r