Solving Schrodinger Equation using time-dependent finite method.

This is a discussion on Solving Schrodinger Equation using time-dependent finite method. within the C Programming forums, part of the General Programming Boards category; Hi there, guys. I have a scholarship at my college, and my teacher asked me to write a program to ...

  1. #1
    Registered User
    Join Date
    Mar 2008
    Posts
    7

    Solving Schrodinger Equation using time-dependent finite method.

    Hi there, guys.

    I have a scholarship at my college, and my teacher asked me to write a program to solve the Schrodinger Equation using the aforementioned method. Mind you, if I don't bring him some results, I may as well lose my scholarship. It's part of a larger project.

    Well, I succeeded in writing a simple program to generate random numbers, write it to a file and make gnuplot read this file and plot those numbers. Now I must incorporate a code that solve the Schrodinger equation, but I have no idea how to do it. The only code (that solves it) that I found so far is in BASIC, and I can't use it as a template because it uses some BASIC-only libraries, and I don't know their counterparts in C.

    Please, some advice?

  2. #2
    Registered User
    Join Date
    Sep 2011
    Posts
    25
    Hey, terminusest.

    Can you link us to the BASIC code for this?

    I'm looking into some articles on solving Schrodinger's equation using the finite method. This is some complex stuff! What is your major? I'm in my third year of Computer / Electrical engineering dual major. Currently I'm reading this article: http://www.physics.udel.edu/~jim/qua...the%20TDSE.pdf It mentions some Euler's method and Fourier analysis. I'm pretty solid on Euler's but we are not going into Fourier until later into this semester. It might help me later on though to get some practice now. Provide a little more info?

  3. #3
    Registered User
    Join Date
    Dec 2012
    Posts
    289
    1. in order for us to help (not do it for you) we need the code you have
    2. BASIC is so "simple" in most of its commands, it is usually pretty easy to change to c coding.

    advice post what you have, so we may further assist you!

  4. #4
    Registered User
    Join Date
    Mar 2008
    Posts
    7
    Quote Originally Posted by Derek Lake View Post
    Hey, terminusest.

    Can you link us to the BASIC code for this?

    I'm looking into some articles on solving Schrodinger's equation using the finite method. This is some complex stuff! What is your major? I'm in my third year of Computer / Electrical engineering dual major. Currently I'm reading this article: http://www.physics.udel.edu/~jim/qua...the%20TDSE.pdf It mentions some Euler's method and Fourier analysis. I'm pretty solid on Euler's but we are not going into Fourier until later into this semester. It might help me later on though to get some practice now. Provide a little more info?
    Here is the BASIC code. My major is Biology, although we in here (Brazil) do not have "majors" and "minors".

    The method used on this code is the split-step (I think), something akin to Euler's. But the method to solve it is not important, maybe the Runge-Kutta 4th order is enough for my purposes.

    Code:
    PROGRAM tdse
    ! motion of a wavepacket incident on a potential
    DIM Re(0:1100),Im(0:1100),Imold(0:1100)
    CALL parameters(x0,k0,width,V0,a,xmin,xmax,n,dx,dt)
    CALL initial_packet(Re(),Im(),x0,k0,width,xmin,n,dx,dt)
    CALL set_up_windows(xmin,xmax,V0,#1,#2,#3,#4)
    CALL draw_potential(V0,a,xmin,xmax,#4)
    DO
    	CALL evolve(Re(),Im(),Imold(),t,V0,a,dx,dx2,dt,xmin,n)
    	CALL plots(Re(),Im(),Imold(),t,xmin,n,dx,#1,#2,#3,#4)
    LOOP until key input
    END
    
    SUB parameters(x0,k0,width,V0,a,xmin,xmax,n,dx,dx2,dt)
    	LET x0 = -15
    	LET width = 1
    	LET k0 = 2
    	LET xmax = 20
    	LET xmin= -xmax
    	LET V0 = 2
    	LET a = 1
    	LET dx = 0.4
    	LET dx2 = dx*dx
    	LET n = (xmax - xmin)/dx
    	LET dt = 0.1
    END SUB
    
    SUB initial_packet(Re(),Im(),x0,k0,width,xmin,n,dx,dt)
    	! initial Gaussian wavepacket
    	LET delta2 = width*width
    	LET A = (2*pi*delta2)^(-0.25)
    	LET b = 0.5*k0*dt
    	FOR i = 1 to n
    		LET x = xmin + (i-1)*dx
    		LET arg = 0.25*(x - x0)^2/delta2
    		LET e = exp(-arg)
    		LET Re(i) = A*cos(k0*(x - x0))*e
    		LET Im(i) = A*sin(k0*(x - x0 - 0.5*b))*e
    	NEXT i
    END SUB
    
    SUB set_up_windows(xmin,xmax,V0,#1,#2,#3,#4)
    	OPEN #1: screen 0,1,0.75,1.0
    	SET WINDOW xmin,xmax,-0.1,0.5
    	PLOT LINES: xmin,0;xmax,0
    	OPEN #2: screen 0,1,0.5,0.75
    	SET WINDOW xmin,xmax,-1,1
    	PLOT LINES: xmin,0;xmax,0
    	OPEN #3: screen 0,1,0.25,0.5
    	SET WINDOW xmin,xmax,-1,1
    	PLOT LINES: xmin,0;xmax,0
    	OPEN #4: screen 0,1,0,0.22
    	SET WINDOW xmin,xmax,-V0,V0
    END SUB
    
    SUB draw_potential(V0,a,xmin,xmax,#4)
    	DECLARE DEF V
    	WINDOW #4
    	SET COLOR "red"
    	PLOT LINES: xmin,0;a,0;a,0;a,V0;a,V0;xmax,V0
    	SET COLOR "black/white"
    	PRINT "total probability = ";
    END SUB
    
    SUB plots(Re(),Im(),Imold(),t,xmin,n,dx,#1,#2,#3,#4)
    	WINDOW #2
    	CLEAR
    	PLOT LINES: xmin,0;xmax,0
    	FOR i = 1 to n
    		LET x = xmin + (i - 1)*dx
    		PLOT x,Im(i);
    	NEXT i
    	WINDOW #1
    	CLEAR
    	PLOT LINES: xmin,0;xmax,0
    	LET Psum = 0
    	FOR i = 1 to n
    		LET x = xmin + (i - 1)*dx
    		LET P = Re(i)*Re(i) + Im(i)*Imold(i)
    		LET Psum = Psum + P*dx
    		PLOT x,P;
    	NEXT i
    	WINDOW #4
    	SET CURSOR 1,20
    	PRINT "        ";
    	SET CURSOR 1,20
    	PRINT Psum
    END SUB
    
    SUB evolve(Re(),Im(),Imold(),t,V0,a,dx,dx2,dt,xmin,n)
      DECLARE DEF V
      FOR i = 1 to n
        LET x = xmin + (i-1)*dx
        LET HIm = V(x,V0,a)*Im(i) -0.5*(Im(i+1) -2*Im(i) +Im(i-1))/dx2
        ! real part defined at multiples of dt
        LET Re(i) = Re(i) + HIm*dt
      NEXT i
      FOR i = 1 to n
        LET x = xmin + (i-1)*dx
        !dt/2 earlier than real part
        LET Imold(i) = Im(i)
        LET HRe = V(x,V0,a)*Re(i) -0.5*(Re(i+1) - 2*Re(i) +Re(i-1)/dx2
        ! dt/2 later than real part
        LET Im(i) = Im(i) - HRe*dt
      NEXT i
      LET t = t + dt 		! time of real part
    END SUB
    
    FUNCTION V(x,V0,a)
    	! step potential
    	IF x > a then
    		LET V = V0
    	ELSE
    		LET V = 0
    	END IF
    END DEF

  5. #5
    Registered User
    Join Date
    Mar 2008
    Posts
    7
    Quote Originally Posted by Crossfire View Post
    1. in order for us to help (not do it for you) we need the code you have
    2. BASIC is so "simple" in most of its commands, it is usually pretty easy to change to c coding.

    advice post what you have, so we may further assist you!
    Here is the little program I wrote that calls a library to speak directly with gnuplot. It had a random number generator for the variables, but I scrapped it. Also, it has some variables that I'm meddling with yet, and serve no purpose for the moment. Sorry for it, I'm a newbie programmer.

    Code:
    #include <stdio.h>
    #include <stdlib.h>
    #include <gnuplot_i.h>
    #include <math.h>
    #include <string.h>
    
    double x = 1.0;
    double y = 7.0;
    double z = 0.0;
    double phase = 0;
    void main(void){
    gnuplot_ctrl *h;
    h = gnuplot_init();
    gnuplot_setstyle(h, "lines");
    FILE *random;
    random = fopen("random.dat","w");
    gnuplot_cmd(h, "set xrange[0:10]; set yrange[0:10]");
      for(x; x < 2; x++){
      fprintf(random,"%d %d\n%d %d\n",rand()%10,rand()%10,rand()%10,rand()%10);
      gnuplot_resetplot(h);
      gnuplot_cmd(h, "plot 'random.dat' using 1:2 with lines");
      x--;
      rewind(random);
      sleep(0.5);
    }
    gnuplot_close(h);
    fclose(random);
    }

  6. #6
    Registered User
    Join Date
    Nov 2012
    Posts
    981
    I was able to run the basic code using TrueBasic free version. It seems like the plotting commands interface to a runtime library only in that version of BASIC.

    To translate the TrueBasic plotting calls to C, you could use a fictional library, "area", which implements the same graph plotting functions as in TrueBasic. A numerical parameter #1 becomes an opaque type AREA*.

    tdse.c:
    Code:
    // PROGRAM tdse
    // motion of a wavepacket incident on a potential
    #include <stdio.h>
    #include <stdlib.h>
    #include <math.h>
    #include "area.h"
    
    const double PI = 3.14159265358979323846264338327950;
    static double Re[1101], Im[1101], Imold[1101];
    void parameters(double *x0, double *k0, double *width, double *V0, double *a, 
            double *xmin, double *xmax, double *n, double *dx, double *dx2, 
            double *dt);
    void initial_packet(double Re[], double Im[], double x0, double k0,
            double width, double xmin, double n, double dx, double dt);
    void set_up_windows(double xmin, double xmax, double V0, 
            AREA **area1, AREA **area2, AREA **area3, AREA **area4);
    void draw_potential(double V0, double a, double xmin, double xmax,
            AREA *area4);
    
    int main(void)
    {
        void evolve(double Re[], double Im[], double Imold[], double *t, 
                double V0, double a, double dx, double dx2, double dt, 
                double xmin, double n);
        void plots(double Re[], double Im[], double Imold[], double *t, 
                double xmin, double n, double dx, AREA *area1, AREA *area2, 
                AREA *area3, AREA *area4);
        double x0, k0, width, V0, a, xmin, xmax, n, dx, dx2, dt;
        parameters(&x0, &k0, &width, &V0, &a, &xmin, &xmax, &n, &dx, &dx2, &dt);
        initial_packet(Re, Im, x0, k0, width, xmin, n, dx, dt);
        AREA *area1, *area2, *area3, *area4;
        set_up_windows(xmin, xmax, V0, &area1, &area2, &area3, &area4);
        draw_potential(V0, a, xmin, xmax, area4);
        double t = 0.0;
        do {
            evolve(Re, Im, Imold, &t, V0, a, dx, dx2, dt, xmin, n);
            plots(Re, Im, Imold, &t, xmin, n, dx, area1, area2, area3, area4);
        } while (getchar() != EOF);
        closearea(area1);
        closearea(area2);
        closearea(area3);
        closearea(area4);    
        return EXIT_SUCCESS;
    }
    
    void parameters(double *x0, double *k0, double *width, double *V0, double *a, 
            double *xmin, double *xmax, double *n, double *dx, double *dx2, 
            double *dt)
    {
        *x0 = -15;
        *width = 1;
        *k0 = 2;
        *xmax = 20;
        *xmin = -*xmax;
        *V0 = 2;
        *a = 1;
        *dx = 0.4;
        *dx2 = (*dx)*(*dx);
        *n = (*xmax - *xmin)/(*dx);
        *dt = 0.1;
    }
    
    void initial_packet(double Re[], double Im[], double x0, double k0,
            double width, double xmin, double n, double dx, double dt)
        // initial Gaussian wavepacket
    {
        double delta2 = width*width;
        double A = pow(2*PI*delta2, -0.25);
        double b = 0.5*k0*dt;
        for (int i=1; i <= n; i++) {
            double x = xmin + (i-1)*dx;
            double arg = 0.25*pow(x - x0, 2)/delta2;
            double e = exp(-arg);
            Re[i] = A*cos(k0*(x - x0))*e;
            Im[i] = A*sin(k0*(x - x0 - 0.5*b))*e;
        }
    }
    
    void set_up_windows(double xmin, double xmax, double V0, 
            AREA **area1, AREA **area2, AREA **area3, AREA **area4)
    {
        *area1 = openarea(0, 1, 0.75, 1.0);
        setwindow(*area1, xmin, xmax, -0.1, 0.5);
        plotline(*area1, xmin, 0, xmax, 0);
        *area2 = openarea(0, 1, 0.5, 0.75);
        setwindow(*area2, xmin, xmax, -1, 1);
        plotline(*area2, xmin, 0, xmax, 0);
        *area3 = openarea(0, 1, 0.25, 0.5);
        setwindow(*area3, xmin, xmax, -1, 1);
        plotline(*area3, xmin, 0, xmax, 0);
        *area4 = openarea(0, 1, 0, 0.22);
        setwindow(*area4, xmin, xmax, -V0, V0);
    }
    
    void draw_potential(double V0, double a, double xmin, double xmax, 
            AREA *area4)
    {
        double V(double x, double V0, double a);
        AREA *area = area4;
        setcolor(area, "red");
        plotline(area, xmin, 0, a, 0);
        plotline(area, a, 0, a, V0);
        plotline(area, a, V0, xmax, V0);
        setcolor(area, "black/white");
        printf("total probability = ");
    }
    
    void plots(double Re[], double Im[], double Imold[], double *t, double xmin, 
            double n, double dx, AREA *area1, AREA *area2, AREA *area3, AREA *area4)
    {
        AREA *area = area2;
        cleararea(area);
        double xmax = 0.0;
        plotline(area, xmin, 0, xmax, 0);
        for (int i=0; i <= n; i++) {
            double x = xmin + (i - 1)*dx;
            plotpoint(area, x, Im[i]);
        }
        area = area1;
        cleararea(area);
        plotline(area, xmin, 0, xmax, 0);
        double Psum = 0;
        for (int i=1 ; i <= n; i++) {
            double x = xmin + (i - 1)*dx;
            double P = Re[i]*Re[i] + Im[i]*Imold[i];
            Psum += P*dx;
            plotpoint(area, x, P);
        }
        area = area4;
        setcursor(area, 1, 20);
        printf("        ");
        setcursor(area, 1, 20);
        printf("%g", Psum);
    }
    
    void evolve(double Re[], double Im[], double Imold[], double *t,
            double V0, double a, double dx, double dx2, double dt,
            double xmin, double n)
    {
        double V(double x, double V0, double a);
        for (int i=1; i <= n; i++) {
            double x = xmin + (i-1)*dx;
            double HIm = V(x,V0,a)*Im[i] - 0.5*(Im[i+1] - 2*Im[i] + Im[i-1]) 
                    / dx2;
            // real part defined at multiples of dt
            Re[i] += HIm*dt;
        }
        for (int i=1; i <= n; i++) {
            double x = xmin + (i-1)*dx;
            // dt/2 earlier than real part
            Imold[i] = Im[i];
            double HRe = V(x,V0,a)*Re[i] - 0.5*(Re[i+1] - 2*Re[i] + Re[i-1])
                    / dx2;
            // dt/2 later than real part
            Im[i] -= HRe*dt;
        }
        *t += dt;        // time of real part
    }
    
    double V(double x, double V0, double a)
    {
        // step potential
        if (x > a) {
            return V0;
        }else{
            return 0;
        }
    }
    If you want the graphics in C, you'll have to implement a compatible set of calls that manipulate the AREA objects. I don't really know how they're supposed to work, so these are just stub functions to get the thing to compile.

    area.h:
    Code:
    #ifndef AREA_H
    #define AREA_H
    
    struct AREA;
    typedef struct AREA AREA;
    
    AREA *openarea(double a, double b, double c, double d);
    void closearea(AREA *area);
    void cleararea(AREA *area);
    void setcolor(AREA *area, const char *colorstr);
    void setcursor(AREA *area, int row, int col);
    void setwindow(AREA *area, double xmin, double xmax, double ymin, 
            double ymax);
    void plotline(AREA *area, double x0, double y0, double x1, double y1);
    void plotpoint(AREA *area, double x0, double y0);
    
    #endif // AREA_H
    area.c:
    Code:
    #include "area.h"
    #include <stdlib.h>
    #include <assert.h>
    
    struct AREA {
        double a;
        double b;
        double c;
        double d;
        void *openfn;
    };
    
    AREA *openarea(double a, double b, double c, double d)
    {
        AREA *area = malloc(sizeof *area);
        assert(area != NULL);
        area->a = a;
        area->b = b;
        area->c = c;
        area->d = d;
        area->openfn = openarea;
        return area;
    }
    
    void closearea(AREA *area)
    {
        assert(area != NULL);
        assert(area->openfn == openarea);
        free(area);
    }
    
    void cleararea(AREA *area) { }
    void setcolor(AREA *area, const char *colorstr) { }
    void setcursor(AREA *area, int row, int col) { }
    void setwindow(AREA *area, double xmin, double xmax, double ymin, double ymax)
    { }
    void plotline(AREA *area, double x0, double y0, double x1, double y1) { }
    void plotpoint(AREA *area, double x0, double y0) { }
    Maybe you replace it with gnuplot calls, but I'm not sure if gnuplot does animation like the TrueBasic runtime calls do.

  7. #7
    TEIAM - problem solved
    Join Date
    Apr 2012
    Location
    Melbourne Australia
    Posts
    1,303
    You shouldn't use assert to test the return value for malloc.
    Code:
    AREA *openarea(double a, double b, double c, double d)
    {
        AREA *area = malloc(sizeof *area);
        assert(area != NULL);
        ...
    The reason is that it is only a debuging tool -> When you are finished debugging, you are supposed to define NDEBUG and all of the asserts are pre-processed out of existence.
    Fact - Beethoven wrote his first symphony in C

  8. #8
    Registered User
    Join Date
    Mar 2008
    Posts
    7
    Quote Originally Posted by c99tutorial View Post
    I was able to run the basic code using TrueBasic free version. It seems like the plotting commands interface to a runtime library only in that version of BASIC.

    To translate the TrueBasic plotting calls to C, you could use a fictional library, "area", which implements the same graph plotting functions as in TrueBasic. A numerical parameter #1 becomes an opaque type AREA*.

    tdse.c:
    Code:
    // PROGRAM tdse
    // motion of a wavepacket incident on a potential
    #include <stdio.h>
    #include <stdlib.h>
    #include <math.h>
    #include "area.h"
    
    const double PI = 3.14159265358979323846264338327950;
    static double Re[1101], Im[1101], Imold[1101];
    void parameters(double *x0, double *k0, double *width, double *V0, double *a, 
            double *xmin, double *xmax, double *n, double *dx, double *dx2, 
            double *dt);
    void initial_packet(double Re[], double Im[], double x0, double k0,
            double width, double xmin, double n, double dx, double dt);
    void set_up_windows(double xmin, double xmax, double V0, 
            AREA **area1, AREA **area2, AREA **area3, AREA **area4);
    void draw_potential(double V0, double a, double xmin, double xmax,
            AREA *area4);
    
    int main(void)
    {
        void evolve(double Re[], double Im[], double Imold[], double *t, 
                double V0, double a, double dx, double dx2, double dt, 
                double xmin, double n);
        void plots(double Re[], double Im[], double Imold[], double *t, 
                double xmin, double n, double dx, AREA *area1, AREA *area2, 
                AREA *area3, AREA *area4);
        double x0, k0, width, V0, a, xmin, xmax, n, dx, dx2, dt;
        parameters(&x0, &k0, &width, &V0, &a, &xmin, &xmax, &n, &dx, &dx2, &dt);
        initial_packet(Re, Im, x0, k0, width, xmin, n, dx, dt);
        AREA *area1, *area2, *area3, *area4;
        set_up_windows(xmin, xmax, V0, &area1, &area2, &area3, &area4);
        draw_potential(V0, a, xmin, xmax, area4);
        double t = 0.0;
        do {
            evolve(Re, Im, Imold, &t, V0, a, dx, dx2, dt, xmin, n);
            plots(Re, Im, Imold, &t, xmin, n, dx, area1, area2, area3, area4);
        } while (getchar() != EOF);
        closearea(area1);
        closearea(area2);
        closearea(area3);
        closearea(area4);    
        return EXIT_SUCCESS;
    }
    
    void parameters(double *x0, double *k0, double *width, double *V0, double *a, 
            double *xmin, double *xmax, double *n, double *dx, double *dx2, 
            double *dt)
    {
        *x0 = -15;
        *width = 1;
        *k0 = 2;
        *xmax = 20;
        *xmin = -*xmax;
        *V0 = 2;
        *a = 1;
        *dx = 0.4;
        *dx2 = (*dx)*(*dx);
        *n = (*xmax - *xmin)/(*dx);
        *dt = 0.1;
    }
    
    void initial_packet(double Re[], double Im[], double x0, double k0,
            double width, double xmin, double n, double dx, double dt)
        // initial Gaussian wavepacket
    {
        double delta2 = width*width;
        double A = pow(2*PI*delta2, -0.25);
        double b = 0.5*k0*dt;
        for (int i=1; i <= n; i++) {
            double x = xmin + (i-1)*dx;
            double arg = 0.25*pow(x - x0, 2)/delta2;
            double e = exp(-arg);
            Re[i] = A*cos(k0*(x - x0))*e;
            Im[i] = A*sin(k0*(x - x0 - 0.5*b))*e;
        }
    }
    
    void set_up_windows(double xmin, double xmax, double V0, 
            AREA **area1, AREA **area2, AREA **area3, AREA **area4)
    {
        *area1 = openarea(0, 1, 0.75, 1.0);
        setwindow(*area1, xmin, xmax, -0.1, 0.5);
        plotline(*area1, xmin, 0, xmax, 0);
        *area2 = openarea(0, 1, 0.5, 0.75);
        setwindow(*area2, xmin, xmax, -1, 1);
        plotline(*area2, xmin, 0, xmax, 0);
        *area3 = openarea(0, 1, 0.25, 0.5);
        setwindow(*area3, xmin, xmax, -1, 1);
        plotline(*area3, xmin, 0, xmax, 0);
        *area4 = openarea(0, 1, 0, 0.22);
        setwindow(*area4, xmin, xmax, -V0, V0);
    }
    
    void draw_potential(double V0, double a, double xmin, double xmax, 
            AREA *area4)
    {
        double V(double x, double V0, double a);
        AREA *area = area4;
        setcolor(area, "red");
        plotline(area, xmin, 0, a, 0);
        plotline(area, a, 0, a, V0);
        plotline(area, a, V0, xmax, V0);
        setcolor(area, "black/white");
        printf("total probability = ");
    }
    
    void plots(double Re[], double Im[], double Imold[], double *t, double xmin, 
            double n, double dx, AREA *area1, AREA *area2, AREA *area3, AREA *area4)
    {
        AREA *area = area2;
        cleararea(area);
        double xmax = 0.0;
        plotline(area, xmin, 0, xmax, 0);
        for (int i=0; i <= n; i++) {
            double x = xmin + (i - 1)*dx;
            plotpoint(area, x, Im[i]);
        }
        area = area1;
        cleararea(area);
        plotline(area, xmin, 0, xmax, 0);
        double Psum = 0;
        for (int i=1 ; i <= n; i++) {
            double x = xmin + (i - 1)*dx;
            double P = Re[i]*Re[i] + Im[i]*Imold[i];
            Psum += P*dx;
            plotpoint(area, x, P);
        }
        area = area4;
        setcursor(area, 1, 20);
        printf("        ");
        setcursor(area, 1, 20);
        printf("%g", Psum);
    }
    
    void evolve(double Re[], double Im[], double Imold[], double *t,
            double V0, double a, double dx, double dx2, double dt,
            double xmin, double n)
    {
        double V(double x, double V0, double a);
        for (int i=1; i <= n; i++) {
            double x = xmin + (i-1)*dx;
            double HIm = V(x,V0,a)*Im[i] - 0.5*(Im[i+1] - 2*Im[i] + Im[i-1]) 
                    / dx2;
            // real part defined at multiples of dt
            Re[i] += HIm*dt;
        }
        for (int i=1; i <= n; i++) {
            double x = xmin + (i-1)*dx;
            // dt/2 earlier than real part
            Imold[i] = Im[i];
            double HRe = V(x,V0,a)*Re[i] - 0.5*(Re[i+1] - 2*Re[i] + Re[i-1])
                    / dx2;
            // dt/2 later than real part
            Im[i] -= HRe*dt;
        }
        *t += dt;        // time of real part
    }
    
    double V(double x, double V0, double a)
    {
        // step potential
        if (x > a) {
            return V0;
        }else{
            return 0;
        }
    }
    If you want the graphics in C, you'll have to implement a compatible set of calls that manipulate the AREA objects. I don't really know how they're supposed to work, so these are just stub functions to get the thing to compile.

    area.h:
    Code:
    #ifndef AREA_H
    #define AREA_H
    
    struct AREA;
    typedef struct AREA AREA;
    
    AREA *openarea(double a, double b, double c, double d);
    void closearea(AREA *area);
    void cleararea(AREA *area);
    void setcolor(AREA *area, const char *colorstr);
    void setcursor(AREA *area, int row, int col);
    void setwindow(AREA *area, double xmin, double xmax, double ymin, 
            double ymax);
    void plotline(AREA *area, double x0, double y0, double x1, double y1);
    void plotpoint(AREA *area, double x0, double y0);
    
    #endif // AREA_H
    area.c:
    Code:
    #include "area.h"
    #include <stdlib.h>
    #include <assert.h>
    
    struct AREA {
        double a;
        double b;
        double c;
        double d;
        void *openfn;
    };
    
    AREA *openarea(double a, double b, double c, double d)
    {
        AREA *area = malloc(sizeof *area);
        assert(area != NULL);
        area->a = a;
        area->b = b;
        area->c = c;
        area->d = d;
        area->openfn = openarea;
        return area;
    }
    
    void closearea(AREA *area)
    {
        assert(area != NULL);
        assert(area->openfn == openarea);
        free(area);
    }
    
    void cleararea(AREA *area) { }
    void setcolor(AREA *area, const char *colorstr) { }
    void setcursor(AREA *area, int row, int col) { }
    void setwindow(AREA *area, double xmin, double xmax, double ymin, double ymax)
    { }
    void plotline(AREA *area, double x0, double y0, double x1, double y1) { }
    void plotpoint(AREA *area, double x0, double y0) { }
    Maybe you replace it with gnuplot calls, but I'm not sure if gnuplot does animation like the TrueBasic runtime calls do.
    Whoa! Calm down! I'm just a newbie, as you can see from my little random number generator program. Ok, let's see if I understand the whole thing: you translated the entire program to C; then, you took the functions and library calls specifics to BASIC and replaced them with a dummy library AREA. I really appreciate your work here, but this is what I asked about in the first place, quoting myself: "The only code (that solves it) that I found so far is in BASIC, and I can't use it as a template because it uses some BASIC-only libraries, and I don't know their counterparts in C."

    Nonetheless, your code is so much cleaner than what I wrote...

    Thank you, not because you solved my former question, but because you were very prestative (don't know if this word exists in english (not sure if there is this word - exists - in english)).

    Then again, some of you guys know which library should I use to replace those BASIC specific calls?

    Also, it would do if I could replace those library calls with printfs to a file, instead of plotting it directly, it would print the resulting numbers in a file, and I could tell gnuplot to use them... I guess...
    Last edited by terminusest; 01-15-2013 at 11:20 AM.

  9. #9
    Registered User
    Join Date
    Nov 2012
    Posts
    981
    Quote Originally Posted by terminusest View Post
    Ok, let's see if I understand the whole thing: you translated the entire program to C; then, you took the functions and library calls specifics to BASIC and replaced them with a dummy library AREA.

    ...

    Nonetheless, your code is so much cleaner than what I wrote...
    It was my intention to make it exactly the same just to bridge the language gap. Didn't you write the BASIC code? Whoever wrote it must know how the method works -- and that somebody is definitely not me. In any case the C version gave me the same numerical answer as the BASIC version.

    Quote Originally Posted by terminusest View Post
    Also, it would do if I could replace those library calls with printfs to a file, instead of plotting it directly, it would print the resulting numbers in a file, and I could tell gnuplot to use them... I guess...
    Yes, then just fill in the AREA functions with printf's that gnuplot will understand it. The only difference I can see is that the BASIC version is using a "clear" function to animate the graph, and I think gnuplot is intended for non-moving plots, but I could be wrong. For example WINDOW #2 followed by CLEAR must involve only the "#2" area, so the C translation becomes
    Code:
     area = area2; // WINDOW #2
     cleararea(area);  // CLEAR
    If you want the graphing done the same way, best would be to sit down with gnuplot and TrueBASIC side by side and look for the appropriate "translations" for each plotting command to accomplish the same behaviour. For example how do you partition of an "area" of the window in gnuplot, how do you "clear" it, how do you draw a line from x0,y0 to x1,y1, etc... if you know these things, then the AREA functions could in principle be used to emulate this functionality of TrueBASIC.

  10. #10
    Registered User
    Join Date
    Mar 2008
    Posts
    7
    The only difference I can see is that the BASIC version is using a "clear" function to animate the graph, and I think gnuplot is intended for non-moving plots, but I could be wrong.
    Well, in fact it does not "animate" in the sense that it refreshes the graph 30 times per second, but as you can see on my little number generator program some posts above that I can clear the graph every half second and plot a new graph.

    If you want the graphing done the same way, best would be to sit down with gnuplot and TrueBASIC side by side and look for the appropriate "translations" for each plotting command to accomplish the same behaviour. For example how do you partition of an "area" of the window in gnuplot, how do you "clear" it, how do you draw a line from x0,y0 to x1,y1, etc... if you know these things, then the AREA functions could in principle be used to emulate this functionality of TrueBASIC.
    I think my limited understanding of programming is what is hindering me from translating the whole thing. I do not fully understand how those functions with multiple arguments - like void parameters(double *x0, double *k0, double *width, double *V0, double *a, double *xmin, double *xmax, double *n, double *dx, double *dx2, double *dt); - really work. I'm studying the code right now, and I'm trying to grasp how can I open a screen on X11 with just C libraries.

    Any clues?

    EDIT: Ok, after some googling I found it's impossible to really "draw" something on console, unless you are after some ASCII art... What about OpenGL, then?

    Advices?
    Last edited by terminusest; 01-16-2013 at 09:53 AM.

  11. #11
    Registered User
    Join Date
    Nov 2012
    Posts
    981
    Your first question seems like its about pointer notation in C. Notice that in the BASIC version you have a function like this

    SUB parameters(x0,k0,width,V0,a,xmin,xmax,n,dx,dx2,dt)
    ...
    END SUB

    Notice that inside the function, the code is allowed to modify the variables. To accomplish the equivalent in C one must use pointers. Basically in practice it just means you prefix everything with * in the parameters list

    void parameters(double *x0, ...);

    And then you must prefix the passed-in local variables with & when you actually call the function.

    parameters(&x0, ...);

    When dealing with pointers many references tend to go into too complicated explanations and start talking about unnecessary things like function pointers and dynamic memory allocation, linked lists, etc., before the basics are covered. I recommend reading this online chapter and trying the examples out to get the hang of it: Practical C Programming, 3rd Edition -- Sample chapter

    Your second question is about drawing primitives. If you're on X11 I suppose you could directly use the X api, but I don't have experience with this and from what I hear it is very painful to use. Gtk has a much better API and is pretty much standard on Linux, and if you include the DLL in your program should work on Windows as well. here is the chapter on drawing primitives: Drawing Primitives

    If you want to use this I would start with a simple example first like make a program which draws a simple plot or something to get the hang of it.

  12. #12
    Registered User
    Join Date
    Mar 2008
    Posts
    7
    Quote Originally Posted by c99tutorial View Post
    Your first question seems like its about pointer notation in C. Notice that in the BASIC version you have a function like this

    SUB parameters(x0,k0,width,V0,a,xmin,xmax,n,dx,dx2,dt)
    ...
    END SUB

    Notice that inside the function, the code is allowed to modify the variables. To accomplish the equivalent in C one must use pointers. Basically in practice it just means you prefix everything with * in the parameters list

    void parameters(double *x0, ...);

    And then you must prefix the passed-in local variables with & when you actually call the function.

    parameters(&x0, ...);

    When dealing with pointers many references tend to go into too complicated explanations and start talking about unnecessary things like function pointers and dynamic memory allocation, linked lists, etc., before the basics are covered. I recommend reading this online chapter and trying the examples out to get the hang of it: Practical C Programming, 3rd Edition -- Sample chapter

    Your second question is about drawing primitives. If you're on X11 I suppose you could directly use the X api, but I don't have experience with this and from what I hear it is very painful to use. Gtk has a much better API and is pretty much standard on Linux, and if you include the DLL in your program should work on Windows as well. here is the chapter on drawing primitives: Drawing Primitives

    If you want to use this I would start with a simple example first like make a program which draws a simple plot or something to get the hang of it.
    I tried GTK before, but couldn't get the thing to compile, GCC always complaining about something or other. On the other hand, OpenGL works here. About pointers, thank you for the materials. I will code something now to see if I can grasp the whole thing. I have some understanding of C, albeit little, but don't know a thing about BASIC.

    Thank you so far, will keep you updated.

  13. #13
    Registered User
    Join Date
    Mar 2008
    Posts
    7
    Well, one weekend after, I could'nt figure out a way. It looks like the code is made to work only with a matching library to draw the graphics. In other words, the whole program is made to use a specific library, and was written for it only.

    I'm trying now to devise the algorithm to solve the equation, and will try to write a new program that uses OpenGL for plotting. Any help would be appreciated. Thank you.

Popular pages Recent additions subscribe to a feed

Similar Threads

  1. Solving a more difficult math equation in c programming
    By karma_jesus45 in forum C Programming
    Replies: 9
    Last Post: 09-14-2011, 06:34 PM
  2. time dependent program
    By shubham in forum C Programming
    Replies: 6
    Last Post: 06-17-2011, 08:31 AM
  3. quadratic equation solving
    By narendrav in forum C Programming
    Replies: 6
    Last Post: 05-18-2011, 08:49 AM
  4. Equation solving
    By Sang-drax in forum A Brief History of Cprogramming.com
    Replies: 3
    Last Post: 11-24-2002, 01:13 PM
  5. equation solving
    By ajlott_coder in forum C Programming
    Replies: 6
    Last Post: 01-26-2002, 01:34 AM

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21