Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Weakly nonlinear theories in electrohydrodynamics

  1. Jun 26, 2010 #1

    hunt_mat

    User Avatar
    Homework Helper

    I derived an equation describing the free surface of an electrified fluid. I am currently seeking travelling wave solutions for this problem, the equation I am looking at is [tex](1-F) f+\frac{1}{90}h^{4}f^{(4)}+\frac{3}{4h}f^{2}-\frac{1}{2}\Bigg( B-\frac{1}{3}\Bigg) hf''+\frac{Ah}{2}\mathcal{H}(f')+\frac{p}{2\rho g}=0[/tex]

    The [tex]\mathcal{H}[/tex] term is the Hilbert transform, [tex]F,h,B,\rho ,g[/tex] are constants

    I have written a numerical program which solves the above equation in MATLAB and it seems to give the sort of solution I am looking for. In order to text this result, if I ignore the nonlinear term, I can derive an analytical solution by the use of Fourier transforms, so if I write:

    [tex]f=\int_{\mathbb{R}}\alpha (k)e^{ikx}dk,\quad p=\int_{\mathbb{R}}\beta (k)e^{ikx}dk[/tex].

    I can get the solution

    [tex]\alpha (k)=-\frac{\beta (k)}{2\rho g}\bigg[ 1-\frac{c}{c_{0}}+\frac{h^{4}k^{4}}{90}+\frac{1}{2}\Bigg( B-\frac{1}{3}\Bigg) h^{2}k^{2}-\frac{Ah|k|}{2}\Bigg]^{-1}[/tex]

    If I compare the above solution to my numerical solution with the nonlinear term ignored, they should, in theory, get the same answer. They don't however. Can anyone see why this should be? The matlab program I use is:

    L=201; %Defines the number of steps, must always be an integer
    x=linspace(-6,6,L); %Defines x
    p=exp(-x.^2); %The pressure distribution
    f=zeros(1,L); %Initial guess for f is f=0;
    Q=0.1; %This the the ratio c/c_0
    %The way forward is to use Newton's method to reduce the equation down to a linear algebra problem and then use Octave
    %to solve the linear algrbra problem.

    %The first step is to compute the Jacobian matrix
    %Do this by diferencing between to close points and dividing through by that difference
    df=10^-10;
    J=zeros(L,L); %This will be the Jacobian.
    it=3; %This will be the number of iterations required to compute the solution
    X=zeros(L,1);
    for k=1:it
    for i=1:L
    for j=1:L
    T=zeros(1,L); %Set all the entries to zero again
    T(i)=df;
    J(i,j)=(g(x,f+T,p,Q,j)-g(x,f,p,Q,j))/df; %Calculates the (i,j)th element of the Jacobian
    end
    end

    for m=1:L
    X(m)=-g(x,f,p,Q,m);
    end

    %Now solve the system of equations
    sol=J\X;
    %sol=inv(J)*X;
    f=f+sol';

    end

    %Now plot the solution!
    plot(x,f);
    xlabel('x');
    ylabel('f(x-ct)');

    The function g is:

    function r=g(x,f,p,Q,i)
    h=1;
    B=0.2;
    E_b=0.1;
    g=9.8065;
    rho=1;
    a_1=1-Q;
    a_2=h^4/90;
    a_3=0.75/h;
    a_4=-0.5*(B-1/3)*h^2;
    a_5=0.5*E_b*h;
    a_6=1/(2*g*rho);
    dx=x(2)-x(1);
    f_ghost=[0 0 f 0 0]; %Add in the ghost cells for f
    N=length(x); %Prepare for the Hilbert transform
    x_half=zeros(1,N-1);
    f_half=zeros(1,N-1);
    for k=1:N-1
    x_half(k)=0.5*(x(k)+x(k+1));
    f_half(k)=0.5*(f(k)+f(k+1));
    end
    grad=gradient(f_half,dx);
    hilbert=trapz(x_half,grad./(x_half-x(i)));
    %hilbert=0; %Start the Hilbert transform at zero
    %for k=2:N-1
    %hilbert=hilbert+(f(k+1)-f(k))/(0.5*(x(k+1)+x(k))-x(i));
    %end
    %Compute the Hilbert transform

    %Now comes the formulae for the finite differencing part of the algorithm.
    b_1=a_1*f_ghost(i+2);
    b_2=(a_2/dx^4)*(f_ghost(i+4)-4*f_ghost(i+3)+6*f_ghost(i+2)-4*f_ghost(i+1)+f_ghost(i)); %Fourth order derivative
    b_3=0; %a_3*f_ghost(i+2)^2; %nonlinear term
    b_4=(a_4/dx^2)*(f_ghost(i+3)-2*f_ghost(i+2)+f_ghost(i+1)); %Second order derivative
    b_5=a_5*hilbert; %Hilbert transform term
    b_6=a_6*p(i); %Pressure term
    r=b_1+b_2+b_3+b_4+b_5+b_6;

    Any idea anyone?
    Mat
     
    Last edited: Jun 26, 2010
  2. jcsd
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook

Can you offer guidance or do you also need help?
Draft saved Draft deleted