#pragma once
class Potential
{
public:
static double Slot1SimpsonsRule(
double a, double A, double B, int n, int nSteps,
double(*fx)(double, double, int));
static double SquareExact(
double a, double b,
double x, double y,
int N);
static void ComputeIntegrals(
double a, int N, int nSteps, std::vector<double>& c);
static double Slot1Exact(
double a, double x, double y,
int N, std::vector<double>& c);
static void SquareFixedPoint(
double a, double b,
int Nx, int Ny, int iterations,
std::vector<std::vector<double>>& U);
static void Slot1FixedPoint(
double a, double b,
int Nx, int Ny, int iterations,
std::vector<std::vector<double>>& U);
static double Disk_f(double theta);
static double Disk_a0(double theta);
static double Disk_an(double c, double theta, int n);
static double Disk_bn(double c, double theta, int n);
static double SimpsonsRule_Disk_a0(int nSteps);
static double SimpsonsRule_Disk_an(
double c, int n, int nSteps);
static double SimpsonsRule_Disk_bn(
double c, int n, int nSteps);
static double Disk_u(
double c, double r, double theta,
double a0, std::vector<double>& an,
std::vector<double>& bn, int Nt,
int nSteps);
};
#include "pch.h"
#include "Potential.h"
double Potential::Slot1SimpsonsRule(
double a, double A, double B, int n, int nSteps,
double(*fx)(double, double, int))
{
double h = (B - A) / nSteps;
double h2 = 2.0 * h;
double s = 0.0;
double t = 0.0;
double x = A + h;
for (int i = 1; i < nSteps; i += 2)
{
s += fx(a, x, n);
x += h2;
}
x = A + h2;
for (int i = 2; i < nSteps; i += 2)
{
t += fx(a, x, n);
x += h2;
}
return h * (fx(a, A, n) + 4 * s + 2 * t + fx(a, B, n)) / 3.0;
}
double Potential::SquareExact(
double a, double b,
double x, double y,
int N)
{
double pi = 4.0 * atan(1.0);
double sum = 0.0;
for (int n = 1; n <= N; n++)
{
double l2 = pi * n / a;
double s1 = sin(0.5 * n * pi) / (n * n);
double s2 = sinh(l2 * y);
double s3 = sinh(l2 * (b - y));
double s4 = sinh(l2 * b);
double s5 = sin(l2 * x);
sum += s1 * ((s2 + s3) / s4) * s5;
}
return 8.0 * sum / (pi * pi);
}
double Slot1Integrand(double a, double x, int n)
{
double pi = 4.0 * atan(1.0);
return exp(-x) * sin(n * pi * x / a);
}
void Potential::ComputeIntegrals(
double a, int N, int nSteps, std::vector<double>& c)
{
c[0] = 0.0;
for (int n = 1; n <= N; n++)
{
c[n] = 2.0 * Slot1SimpsonsRule(a, 0.0, a,
n, nSteps, Slot1Integrand) / a;
}
}
double Potential::Slot1Exact(
double a, double x, double y,
int N, std::vector<double>& c)
{
double pi = 4.0 * atan(1.0), uxy = 0.0;
for (int n = 1; n <= N; n++)
uxy += c[n] * sin(n * pi * x / a) * exp(-n * pi * y / a);
return uxy;
}
static double SquareF(double x)
{
if (x > 0.0 && x < 0.5)
return 2.0 * x;
else if (x >= 0.5 && x < 1.0)
return 2.0 - 2.0 * x;
else
return 0.0;
}
void Potential::SquareFixedPoint(
double a, double b,
int Nx, int Ny, int iterations,
std::vector<std::vector<double>>& U)
{
double deltaX = a / Nx, deltaY = b / Ny;
U.resize(Nx + 1LL);
for (int i = 0; i <= Nx; i++)
U[i].resize(Ny + 1LL);
for (int i = 0; i <= Ny; i++)
U[0][i] = 0.0;
for (int i = 1; i <= Ny; i++)
U[Nx][i] = 0.0;
for (int i = 0; i <= Nx; i++)
U[i][0] = SquareF(i * deltaX);
for (int i = 0; i <= Nx; i++)
U[i][Ny] = SquareF(i * deltaX);
for (int m = 0; m < iterations; m++)
{
for (int j = 1; j < Nx; j++)
{
for (int k = 1; k < Ny; k++)
U[j][k] = 0.25 *
(U[j + 1][k] + U[j - 1][k] +
U[j][k + 1] + U[j][k - 1]);
}
}
}
void Potential::Slot1FixedPoint(
double a, double b,
int Nx, int Ny, int iterations,
std::vector<std::vector<double>>& U)
{
double deltaX = a / Nx, deltaY = b / Ny;
U.resize(Nx + 1);
for (int i = 0; i <= Nx; i++)
U[i].resize(Ny + 1, 0.0);
for (int i = 0; i <= Nx; i++)
U[i][0] = exp(-i * deltaX);
for (int i = 0; i < iterations; i++)
{
for (int j = 1; j < Nx; j++)
{
for (int k = 1; k < Ny; k++)
U[j][k] = 0.25 *
(U[j + 1][k] + U[j - 1][k] +
U[j][k + 1] + U[j][k - 1]);
}
}
}
double Potential::Disk_f(double theta)
{
return theta * theta;
}
double Potential::Disk_a0(double theta)
{
double pi = 4.0 * atan(1.0);
return Disk_f(theta);
}
double Potential::Disk_an(double c, double theta, int n)
{
double pi = 4.0 * atan(1.0);
return Disk_f(theta) * cos(n * theta) / (pi * pow(c, n));
}
double Potential::Disk_bn(double c, double theta, int n)
{
double pi = 4.0 * atan(1.0);
return Disk_f(theta) * cos(n * theta) / (pi * pow(c, n));
}
double Potential::SimpsonsRule_Disk_a0(int nSteps)
{
double pi = 4.0 * atan(1.0);
double h = (pi + pi) / nSteps;
double h2 = 2.0 * h;
double s = 0.0;
double t = 0.0;
double theta = -pi + h;
for (int i = 1; i < nSteps; i += 2)
{
s += Disk_a0(theta);
theta += h2;
}
theta = -pi + h2;
for (int i = 2; i < nSteps; i += 2)
{
t += Disk_a0(theta);
theta += h2;
}
return h * (Disk_a0(0.0) + 4 * s + 2 * t + Disk_a0(pi)) / 3.0;
}
double Potential::SimpsonsRule_Disk_an(
double c, int n, int nSteps)
{
double pi = 4.0 * atan(1.0);
double h = (pi + pi) / nSteps;
double h2 = 2.0 * h;
double s = 0.0;
double t = 0.0;
double theta = -pi + h;
for (int i = 1; i < nSteps; i += 2)
{
s += Disk_an(c, theta, n);
theta += h2;
}
theta = -pi + h2;
for (int i = 2; i < nSteps; i += 2)
{
t += Disk_an(c, theta, n);
theta += h2;
}
double integral =
h * (Disk_an(c, 0.0, n) + 4 * s + 2 * t +
Disk_an(c, pi, n)) / 3.0;
return integral / (pi * pow(c, n));
}
double Potential::SimpsonsRule_Disk_bn(
double c, int n, int nSteps)
{
double pi = 4.0 * atan(1.0);
double h = (pi + pi) / nSteps;
double h2 = 2.0 * h;
double s = 0.0;
double t = 0.0;
double theta = -pi + h;
for (int i = 1; i < nSteps; i += 2)
{
s += Disk_bn(c, theta, n);
theta += h2;
}
theta = -pi + h2;
for (int i = 2; i < nSteps; i += 2)
{
t += Disk_bn(c, theta, n);
theta += h2;
}
double integral =
h * (Disk_bn(c, 0.0, n) + 4 * s + 2 * t +
Disk_bn(c, pi, n)) / 3.0;
return integral / (pi * pow(c, n));
}
double Potential::Disk_u(
double c, double r, double theta,
double a0, std::vector<double>& an,
std::vector<double>& bn, int Nt,
int nSteps)
{
double urt = a0;
for (int n = 1; n <= Nt; n++)
{
urt += pow(r, n) * (
an[n] * cos(n * theta) +
bn[n] * sin(n * theta));
}
return urt;
}
Author: jamespatewilliamsjr
My whole legal name is James Pate Williams, Jr. I was born in LaGrange, Georgia approximately 70 years ago. I barely graduated from LaGrange High School with low marks in June 1971. Later in June 1979, I graduated from LaGrange College with a Bachelor of Arts in Chemistry with a little over a 3 out 4 Grade Point Average (GPA). In the Spring Quarter of 1978, I taught myself how to program a Texas Instruments desktop programmable calculator and in the Summer Quarter of 1978 I taught myself Dayton BASIC (Beginner's All-purpose Symbolic Instruction Code) on LaGrange College's Data General Eclipse minicomputer. I took courses in BASIC in the Fall Quarter of 1978 and FORTRAN IV (Formula Translator IV) in the Winter Quarter of 1979. Professor Kenneth Cooper, a genius poly-scientist taught me a course in the Intel 8085 microprocessor architecture and assembly and machine language. We would hand assemble our programs and insert the resulting machine code into our crude wooden box computer which was designed and built by Professor Cooper. From 1990 to 1994 I earned a Bachelor of Science in Computer Science from LaGrange College. I had a 4 out of 4 GPA in the period 1990 to 1994. I took courses in C, COBOL, and Pascal during my BS work. After graduating from LaGrange College a second time in May 1994, I taught myself C++. In December 1995, I started using the Internet and taught myself client-server programming. I created a website in 1997 which had C and C# implementations of algorithms from the "Handbook of Applied Cryptography" by Alfred J. Menezes, et. al., and some other cryptography and number theory textbooks and treatises.
View all posts by jamespatewilliamsjr
