CAnalyticDike.cpp

Go to the documentation of this file.

#include "CAnalyticDike.h"
#include <math.h>
#include <iostream>

CAnalyticDike::CAnalyticDike()
{}
CAnalyticDike::~CAnalyticDike()
{}

void CAnalyticDike::CalcResponse(double freq, double distance)
{
        double sigma1 = 1./(rho1*1e11);
        double sigma2 = 1./(rho2*1e11);
        double omega = 2 * PI * freq;
        double x = distance * 100;
        double l = length * 100;
        double h = height * 100;
        dcomp U11, U21;
        dcomp U1n, U2n;
        dcomp U1,U2;
        dcomp k1,k2;
        double d1 = sqrt(1./(2 * PI * sigma1 * omega));
        double d2 = sqrt(1./(2 * PI * sigma2 * omega));
        double p1 = d1/ sqrt(2.);
        double p2 = d2/sqrt(2.);
        double h1 = h / p1;
        double h2 = h / p2;
        int nmax = 20;
        cout << "p1: " << p1 << " p2 " << p2 << " h1: " << h1 << " h2 " << h2 << endl;
        if ( sub ==0)
        {
                U1 = 0;
                U2 = 0;
                for (int n = 1; n < nmax; ++n)
                {
                        k1 = sqrt(pow(n*PI,2) + I * pow(h1,2));
                        k2 = sqrt(pow(n*PI,2) + I * pow(h2,2)); 
                        U1n = pow(1.*n,2) * cosh( k1 * x/h) / (pow(k1*k2,2) * (cosh(k1 * l/(2*h)) 
                                + k1/k2 * pow(h2/h1,2)*sinh(k1 * l/(2*h)))); 
                        U2n = pow(1.*n,2) * exp(k2 *l/(2*h)) * sinh(k1 * l/(2*h)) * exp(- k2 * x/(2*h)) /
                        (pow(k1*k2,2)*(sinh(k1*l/(2*h))+ k2/k1 * pow(h1/h2,2)*cosh(k1*l/(2*h))));
                        //U2n = pow(1.*n,2) * exp(k2*(l-2*x)/(2*h)) * sinh(k1 * l/(2*h)) /
                        //(pow(k1*k2,2)*(sinh(k1*l/(2*h))+ k2/k1 * pow(h1/h2,2)*cosh(k1*l/(2*h))));
                        U1 += U1n;
                        U2 += U2n;
                        cout << "U1diff: " << U1n/U1 << " U2diff: " << U2n/U2 << endl;
                }       
                //cout << "k1: " << k1 << " k2: " << k2 << " U1: " << U1 << " U2 " << U2 << endl;
                Z1 = omega * h/h1 * ( sqrt(I) / tanh(sqrt(I)*h1) - 2 * pow(PI,2) * I/h1 * 
                        (pow(h1,2) - pow(h2,2)) * U1); 
                Z2 = omega * h/h2 * ( sqrt(I) / tanh(sqrt(I)*h2) - 2 * pow(PI,2) * I/h2 * 
                        (pow(h2,2) - pow(h1,2)) * U2); 
        }
        else
        {
        }
}

---