diff --git a/Hausaufgaben/ComputerphysikAbgabe/src/ComputerphysikAbgabe.c b/Hausaufgaben/ComputerphysikAbgabe/src/ComputerphysikAbgabe.c index 33b8a42..99751fe 100644 --- a/Hausaufgaben/ComputerphysikAbgabe/src/ComputerphysikAbgabe.c +++ b/Hausaufgaben/ComputerphysikAbgabe/src/ComputerphysikAbgabe.c @@ -1,15 +1,18 @@ #include #include -#include +#include int mmax = 30; //globale Variable, Maximale Anzahl der Schritte -double InfInt(int n, double a, double z,double eps,double (*func)(double, double, double)); //InfInt wird bei Potenzial gebraucht deshalb ist es hier definiert -double f1(double x, double a, double z) { // Zu integrierende Funktion + +double f1(double x, double a, double z) // Zu integrierende Funktion +{ return (exp(-(x * x) / (a * a)) / sqrt(x * x + z * z)); } -double f2(double x, double a, double z){ - return((exp(-(x*x)/(a*a))*z)/pow((x*x)+(z*z),1.5)); + +double f2(double x, double a, double z) { + return ((exp(-(x * x) / (a * a)) * z) / pow((x * x) + (z * z), 1.5)); } + double trapez(int n, double a, double z, double (*func)(double, double, double), double x1, double x2) { //Trapeze double sum; @@ -34,11 +37,8 @@ double IntRomb(double x1, double x2, int n, double a, double z, for (i = 0; i <= mmax; i++) { T[i] = (double*) malloc(((mmax + 1)) * sizeof(double)); //Speicherplatz für 2-dim Array allozieren - } - h[0] = (x2 - x1) / (n); //Erste Schrittweite berechnen - T[0][0] = trapez(n, a, z, func, x1, x2); //erstes Trapez berechnen for (i = 1; i <= mmax; i++) { @@ -52,17 +52,14 @@ double IntRomb(double x1, double x2, int n, double a, double z, / ((h[i] * h[i]) - (h[i - k] * h[i - k])) * T[i][k - 1]) + (h[i] * h[i] / (h[i] * h[i] - h[i - k] * h[i - k]) * T[i - 1][k - 1]); - } - //printf("%d",i); if (abs(T[i][i] - T[i - 1][i - 1]) <= eps) { //Abschätzen ob das Ergebnis genau genug ist wenn ja aus Schleife raus ergebnis = T[i][i]; break; } - } - if (i==mmax){ + if (i == mmax) { printf("Fehler, Integral ist ungenau"); // Fehler falls Abbruchbedingung der vorherigen Schleife nicht greift } //Speicherplatz befreien @@ -74,85 +71,134 @@ double IntRomb(double x1, double x2, int n, double a, double z, return ergebnis; //Rückgabewert ist Ergebnis } -double potenzial(int n, double a, double z, - double (*func)(double, double, double), double eps) { - double V; //Variable Potential - double alpha = 1/sqrt(M_PI); //normierung alpha - double Q = 4; //Gesamtladung Q - V = alpha * (Q / a) * InfInt(n, a, z, eps,func); //Potential berechnen - return V; - -} -double InfInt(int n, double a, double z, double eps,double (*func)(double, double, double)) { +double InfInt(int n, double a, double z, double eps, + double (*func)(double, double, double)) { double x1 = -1; //Grenzen definieren - double x2 = 0; - double R=IntRomb(x1,x2,n,a,z,func,eps); //Variable zum speichern des Ergebnis, Integral für erste Grenzen ausrechnen + double x2 = 0; //hier ist x2 null, da das Integral achsensymmetrisch ist + double R = IntRomb(x1, x2, n, a, z, func, eps); //Variable zum speichern des Ergebnis, Integral für erste Grenzen ausrechnen double lastR; //Vergleichsvariable definieren do { - lastR = R;//letztes Ergebnis speichern - x1 = x1*2; //Grenzen vergrößern - x2 = x2*2; + lastR = R; //letztes Ergebnis speichern + x1 = x1 * 2; //Grenzen vergrößern + x2 = x2 * 2; R = IntRomb(x1, x2, n, a, z, func, eps); //neues Ergebnis berechnen - } while ((lastR/R<(1-eps)) || (lastR/R>(1+eps))); //Integralwerte vergleichen ,sodass deren Änderung kleiner als eps ist - return 2*R; + } while ((lastR / R < (1 - eps)) || (lastR / R > (1 + eps))); //Integralwerte vergleichen ,sodass deren Änderung kleiner als eps ist + return 2 * R; //auf Grund von achsensymmetrie 2 mal } -void efeld(double a, double z, int n, double eps, - double *dfp,double (*func)(double, double, double)) { + +double potenzial(int n, double a, double z, double Q, + double (*func)(double, double, double), double eps) { + double V; //Variable Potential + double alpha = 1 / sqrt(3.141); //normierung alpha + V = alpha * (Q / a) * InfInt(n, a, z, eps, func); //Potential berechnen + return V; +} + +double ges_potenzial(int n, double a1, double a2, double z, double Q1, + double Q2, double d, double (*func)(double, double, double), double eps) { /*Berechnung des gesamten Potenzials der beiden Drähte durch Superposition, wobei für die + z-Koordinate, bzw den Abstand der Ladung zum zweiten Drahtes z-d gilt.*/ + double V_ges; + V_ges = potenzial(n, a1, z, Q1, f1, eps) + + potenzial(n, a2, z - d, Q2, f1, eps); + return V_ges; +} + +void efeld(double a, double z, int n, double Q, double eps, double *dfp, + double (*func)(double, double, double)) { double f[2]; //Array zum speichern der Funktionswerte - double h =0.5* z; //h relativ groß wird in Schleife verkleinert + double h = 0.5 * z; //h relativ groß wird in Schleife verkleinert double lastdf; - f[0] = potenzial( n, a, z +(h), func, eps); //f(z+h) berechnen - f[1] = potenzial(n, a, z - h, func, eps); //f(z-h) berechnen - *dfp = -1* (f[0] - f[1]) / (double)(2 * h); + f[0] = potenzial(n, a, z + (h), Q, func, eps); //f(z+h) berechnen + f[1] = potenzial(n, a, z - h, Q, func, eps); //f(z-h) berechnen + *dfp = -1 * (f[0] - f[1]) / (double) (2 * h); do { - h=h/2; //h kleiner machen - lastdf=*dfp; //Alte Ableitung auf neuer Variable speichern - f[0] = potenzial( n, a, z +(h), func, eps); //f(z+h) berechnen - f[1] = potenzial(n, a, z - h, func, eps); //f(z-h) berechnen + h = h / 2; //h kleiner machen + lastdf = *dfp; //Alte Ableitung auf neuer Variable speichern + f[0] = potenzial(n, a, z + (h), Q, func, eps); //f(z+h) berechnen + f[1] = potenzial(n, a, z - h, Q, func, eps); //f(z-h) berechnen + *dfp = -1 * (f[0] - f[1]) / (2 * h); //Differenzenquotient + } while ((lastdf / (*dfp) < (1 - eps)) || (lastdf / (*dfp) > (1 + eps))); //Werte vergleichen ,sodass deren Änderung kleiner als eps ist +} - *dfp = -1* (f[0] - f[1]) / (2 * h); //Differenzenquotient - } while((lastdf/(*dfp)<(1-eps)) || (lastdf/(*dfp)>(1+eps))); //Werte vergleichen ,sodass deren Änderung kleiner als eps ist +double ges_efeld(int n, double a1, double a2, double z, double Q1, double Q2, + double d, double *df1, double *df2, + double (*func)(double, double, double), double eps) { //verwendet numerisch berechnetes E-Feld + efeld(a1, z, n, Q1, eps, df1, func); //Berechnung E-Feld des 1. Drahtes + efeld(a2, d - z, n, Q2, eps, df2, func); //Berechnung E-Feld des 2. Drahtes + double E_ges = *df1 + *df2; //Superposition der E-Felder der einzelnen Drähte + return E_ges; +} +double secant_pauline(int n, double a1, double a2, double Q1, double Q2, + double d, double *df1, double *df2, + double (*func)(double, double, double), double eps) { + //Setze die startgrenzen ein bisschen von den Drähten entfernt, da das Feld auf den Drähten null ist + double g1 = eps; + double g2 = d - eps; + + double zero; + + //Calculating secants until the zero changes less than eps + while ((g1 / g2) > (1 + eps) || (g1 / g2) < (1 - eps)) { + //Berechnen von mx + c aus Grenzen p1 und p2 + double m = (ges_efeld(n, a1, a2, g2, Q1, Q2, d, df1, df2, func, eps) + - ges_efeld(n, a1, a2, g1, Q1, Q2, d, df1, df2, func, eps)) + / (g2 - g1); + double c = ges_efeld(n, a1, a2, g1, Q1, Q2, d, df1, df2, func, eps) + - (g1 * m); + //Mit dem Sekantenverfahren die nullstelle berechnen + zero = (-1) * (c / m); + g2 = g1; + g1 = zero; + } + //returning the result + return zero; } + int main(void) { - double z, a; - int n=30; - double eps = 1e-8; //genauigkeit epsilon + double z, a1, a2; + int n = 30; + const double eps = 1e-8; //genauigkeit epsilon + const double eps2 = 1e-4; + double d = 12; //Abstand 2.Draht zum 1. + double Q1 = 1; //Gesamtladung der Drähte 1 und 2 + double Q2 = 4; + double df; + double df1; //Variable zum speichern des Differenzenquotienten + double df2; - double df; //Variable zum speichern des Differenzenquotienten //Eingabe der Variablen: - printf("Bitte geben Sie den Abstand der Punktladung von der x-Achse z ein\n "); + printf( + "Bitte geben Sie den Abstand der Punktladung von der x-Achse z ein\n "); scanf("%lf", &z); - printf("Bitte geben Sie a ein\n"); - scanf("%lf", &a); - - - - printf("Unendliches Inegral %lf\n", InfInt(n, a, z, eps,f1)); //Integral berechnen mit InfInt - printf("Das Potenzial ist %f\n", potenzial(n,a,z,f1,eps)); //Potenzial aus Integral berechnen + printf("Bitte geben Sie a1 ein (0 < a1 <= 4):\n"); //Eingabe der Breiten der Drähte + scanf("%lf", &a1); + printf("Bitte geben Sie a2 ein (a1/10 <= a2 <= a1):\n"); + scanf("%lf", &a2); + //H1 berechnung des Potentials + printf("Das Potenzial ist %f\n", potenzial(n, a1, z, Q1, f1, eps)); //Potenzial aus Integral berechnen //H2.1 Berechnung des Efeldes mit Hilfe numerischer Ableitung - efeld(a,z,n,eps,&df,f1); //E-feld berechnen + efeld(a1, z, n, Q1, eps, &df, f1); //E-feld berechnen printf("Das E-Feld an der Stelle %3.2lf ist:%f\n", z, df); //H2.2 Berechnung des Efeldes mit Hilfe analytischer Ableitung - printf("Das E-Feld an der Stelle %3.2lf ist:%f\n", z, potenzial(n,a,z,f2,eps)); - double diff=fabs(df - potenzial(n,a,z,f2,eps)); - printf("Differenz %1.10f", diff); - printf("Achtung schleife"); - for(int i=0; i<20; i++){ - efeld(a,z,n,eps,&df,f1); - printf("Das E-Feld an der Stelle %3.2lf ist:%f\n", z, df); - printf("Das E-Feld an der Stelle %3.2lf ist:%f\n", z, potenzial(n,a,z,f2,eps)); - double diff=fabs(df - potenzial(n,a,z,f2,eps)); - printf("ZZ %3.2lf\n", z); - printf("Differenz %1.10f \n", diff); - - z=z*2; - - } + printf("Das E-Feld an der Stelle %3.2lf ist:%f\n", z, + potenzial(n, a1, z, Q1, f2, eps)); + double diff = fabs(df - potenzial(n, a1, z, Q1, f2, eps)); + printf("Differenz %1.10f", diff); + //H3 zwei Drähte: Setzte a1 als 4 + a1 = 4; + printf("Gesamtes Potenzial der beiden Drähte: %f\n", + ges_potenzial(n, a1, a2, z, Q1, Q2, d, f1, eps2)); + printf("Gesamtes E-Feld (numerisch berechnet) durch Addition ist: %f\n", + ges_efeld(n, a1, a2, z, Q1, Q2, d, &df1, &df2, f1, eps2)); + + printf("Nullstelle %f\n", + secant_pauline(n, a1, a2, Q1, Q2, d, &df1, &df2, f1, eps2)); + z = 0.1; }