/* osci.cpp: Solution of the one-dimensional Schrodinger equation fora particle in a harmonic potential, using the shooting method.To compile and link with gnu compiler, type: g++ -o osci osci.cppTo run the current C++ program, simply type: osciPlot by gnuplot: /GNUPLOT> set terminal windows/GNUPLOT> plot "psi-osc.dat" with lines */#include <stdio.h>#include <stdlib.h>#include <math.h>#define MAX(a, b) (((a) > (b)) ? (a) : (b))int main(int argc, char*argv[]){// Runtime constantsconst static double Epsilon = 1e-10; // Defines the precision of//... energy calculations//APPENDIX C. C++ CODE FOR SOLVING THE SCHRODINGER EQUATIONconst static int N_of_Divisions = 1000;const static int N_max = 5; //Number of calculated EigenstatesFILE *Wavefunction_file, *Energy_file, *Potential_file;Wavefunction_file = fopen("psi-osc.dat","w");Energy_file = fopen("E_n_Oszillator.dat","w");Potential_file = fopen("HarmonicPotentialNoDim.dat", "w");if (!(Wavefunction_file && Energy_file && Potential_file)){ printf("Problems to create files output.\n"); exit(2); }/* Physical parameters using dimensionless quantities.ATTENTION: We set initially: hbar = m = omega = a = 1, andreintroduce physical values at the end. According to Eq.(4.117),the ground state energy then is E_n = 0.5. Since the wave functionvanishes only at -infinity and +infinity, we have to cut off thecalculation somewhere, as given by ’xRange’. If xRange is chosentoo large, the open (positive) end of the wave function candiverge numerically in this simple shooting approach. */const static double xRange = 12; // xRange=11.834 corresponds to a//... physical range of -20fm < x < +20fm, see after Eq.(4.199).const static double h_0 = xRange / N_of_Divisions;double* E_pot = new double[N_of_Divisions+1];double dist;for (int i = 0; i <= N_of_Divisions; ++i){ // Harmonic potential, as given in Eq. (4.115), but dimensionlessdist = i*h_0 - 0.5*xRange;E_pot[i] = 0.5*dist*dist; // E_pot[i]=0;//E_pot=0:Infinite Well!fprintf(Potential_file, "%16.12e \t\t %16.12e\n", dist, E_pot[i]);}fclose(Potential_file);/* Since the Schrodinger equation is linear, the amplitude of thewavefunction will be fixed by normalization.At left we set it small but nonzero. */const static double Psi_left = 1.0e-3; // left boundary conditionconst static double Psi_right = 0.0; // right boundary conditiondouble *Psi, *EigenEnergies;// Arrays to hold the resultsPsi = new double[N_of_Divisions+1]; //N_of_Points = N_of_Divisions+1EigenEnergies = new double[N_max+1];Psi[0] = Psi_left;Psi[1] = Psi_left + 1.0e-3; // Add arbitrary small valueint N_quantum;//N_quantum is Energy Quantum Numberint Nodes_plus; // Number of nodes (+1) in wavefunctiondouble K_square;// Square of wave vector// Initial Eigen-energy search limitsdouble E_lowerLimit = 0.0;// Eigen-energy must be positivedouble E_upperLimit = 10.0;int End_sign = -1;char Limits_are_defined ="false";double Normalization_coefficient;double E_trial;// MAIN LOOP begins:-----------------------------------for(N_quantum=1; N_quantum <= N_max; ++N_quantum){// Find the eigen-values for energy. See theorems (4.1) and (4.2).Limits_are_defined ="false";while (Limits_are_defined == false){ /* First, determine an upper limit for energy, so that the wavefunctionPsi[i] has one node more than physically needed. */Nodes_plus = 0;E_trial = E_upperLimit;for (int i=2; i <= N_of_Divisions; ++i){ K_square = 2.0*(E_trial - E_pot[i]);// Now use the NUMEROV-equation (4.197) to calculate wavefunctionPsi[i] = 2.0*Psi[i-1]*(1.0 - (5.0*h_0*h_0*K_square / 12.0))/(1.0 + (h_0*h_0*K_square/12.0))-Psi[i-2];if (Psi[i]*Psi[i-1] < 0) ++Nodes_plus;}/* If one runs into the following condition, the modificationof the upper limit was too aggressive. */if (E_upperLimit < E_lowerLimit)E_upperLimit = MAX(2*E_upperLimit, -2*E_upperLimit);if (Nodes_plus > N_quantum) E_upperLimit *= 0.7;else if (Nodes_plus < N_quantum) E_upperLimit *= 2.0;else Limits_are_defined =" true"; // At least one node should appear.} // End of the loop: while (Limits_are_defined == false)// Refine the energy by satisfying the right boundary condition.End_sign = -End_sign;while ((E_upperLimit - E_lowerLimit) > Epsilon){ E_trial = (E_upperLimit + E_lowerLimit) / 2.0;for (int i=2; i <= N_of_Divisions; ++i){ // Again eq.(4.197) of the Numerov-algorithm:K_square = 2.0*(E_trial - E_pot[i]);Psi[i] = 2.0*Psi[i-1] * (1.0 - (5.0*h_0*h_0*K_square / 12.0))/(1.0 + (h_0*h_0*K_square/12.0))-Psi[i-2];}if (End_sign*Psi[N_of_Divisions] > Psi_right) E_lowerLimit = E_trial;else E_upperLimit = E_trial;} // End of loop: while ((E_upperLimit - E_lowerLimit) > Epsilon)//APPENDIX C. C++ CODE FOR SOLVING THE SCHRÖDINGER EQUATION// Initialization for the next iteration in main loopE_trial = (E_upperLimit+E_lowerLimit)/2;EigenEnergies[N_quantum] = E_trial;E_upperLimit = E_trial;E_lowerLimit = E_trial;// Now find the normalization coefficientdouble Integral = 0.0;for (int i=1; i <= N_of_Divisions; ++i){ // Simple integrationIntegral += 0.5*h_0*(Psi[i-1]*Psi[i-1]+Psi[i]*Psi[i]);}Normalization_coefficient = sqrt(1.0/Integral);// Output of normalized dimensionless wave functionfor (int i=0; i <=N_of_Divisions; ++i){ fprintf(Wavefunction_file, "%16.12e \t\t %16.12e\n",i*h_0 - 0.5*xRange, Normalization_coefficient*Psi[i]);}fprintf(Wavefunction_file,"\n");} // End of MAIN LOOP. --------------------------------fclose(Wavefunction_file);/*Finally convert dimensionless units in real units. Note thatenergy does not depend explicitly on the particle’s mass anymore:hbar = 1.05457e-34;// Planck constant/2piomega = 5.34e21; // Frequency in 1/sMeV = 1.602176487e-13; // in JThe correct normalization would be hbar*omega/MeV = 3.5148461144,but we use the approximation 3.5 for energy-scale as in chap. 4.9 */const static double Energyscale = 3.5;// in MeV// Output with rescaled dimensions; assign Energy_fileprintf("Quantum Harmonic Oscillator, program osci.cpp\n");printf("Energies in MeV:\n");printf("n \t\t E_n\n");for (N_quantum=1; N_quantum <= N_max; ++N_quantum){ fprintf(Energy_file,"%d \t\t %16.12e\n", N_quantum-1,Energyscale*EigenEnergies[N_quantum]);printf("%d \t\t %16.12e\n", N_quantum-1,Energyscale*EigenEnergies[N_quantum]);}fprintf(Energy_file,"\n");fclose(Energy_file);printf("Wave-Functions in File: psi_osc.dat \n");printf("\n");return 0;}

I tried to compile the code and found only one compiler error:(CODE)This line tries to initialize a single char with a string what's not possible. You should change this line to use a bool and remove

This line tries to initialize a single char with a string what's not possible. You should change this line to use a bool and remove "-s around all occurances of true and false in the rest of the code:

After changing this I still got warnings/errors, but they only are about using function fopen which is considered to be unsafe, you should use fopen_s instead.

Hope that helps,

ZOPPO

PS: I agree with phoffric, you should at least post the error message with this kind of questions, and indenting the code (and even better using a CODE-blocks) makes it much easier to read/understand such a question.

i used my IDE to indent your code and did what Andy has suggested:

#define _CRT_SECURE_NO_WARNINGS // only needed for VC compiler/* osci.cpp: Solution of the one-dimensional Schrodinger equation fora particle in a harmonic potential, using the shooting method.To compile and link with gnu compiler, type: g++ -o osci osci.cppTo run the current C++ program, simply type: osciPlot by gnuplot: /GNUPLOT> set terminal windows/GNUPLOT> plot "psi-osc.dat" with lines */#include <stdio.h>#include <stdlib.h>#include <math.h>#define MAX(a, b) (((a) > (b)) ? (a) : (b))int main(int argc, char*argv[]){// Runtime constants const static double Epsilon = 1e-10; // Defines the precision of //... energy calculations //APPENDIX C. C++ CODE FOR SOLVING THE SCHRODINGER EQUATION const static int N_of_Divisions = 1000; const static int N_max = 5; //Number of calculated Eigenstates FILE *Wavefunction_file, *Energy_file, *Potential_file; Wavefunction_file = fopen("psi-osc.dat","w"); Energy_file = fopen("E_n_Oszillator.dat","w"); Potential_file = fopen("HarmonicPotentialNoDim.dat", "w"); if (!(Wavefunction_file && Energy_file && Potential_file)) { printf("Problems to create files output.\n"); exit(2); } /* Physical parameters using dimensionless quantities. ATTENTION: We set initially: hbar = m = omega = a = 1, and reintroduce physical values at the end. According to Eq.(4.117), the ground state energy then is E_n = 0.5. Since the wave function vanishes only at -infinity and +infinity, we have to cut off the calculation somewhere, as given by ’xRange’. If xRange is chosen too large, the open (positive) end of the wave function can diverge numerically in this simple shooting approach. */ const static double xRange = 12; // xRange=11.834 corresponds to a //... physical range of -20fm < x < +20fm, see after Eq.(4.199). const static double h_0 = xRange / N_of_Divisions; double* E_pot = new double[N_of_Divisions+1]; double dist; for (int i = 0; i <= N_of_Divisions; ++i) { // Harmonic potential, as given in Eq. (4.115), but dimensionless dist = i*h_0 - 0.5*xRange; E_pot[i] = 0.5*dist*dist; // E_pot[i]=0;//E_pot=0:Infinite Well! fprintf(Potential_file, "%16.12e \t\t %16.12e\n", dist, E_pot[i]); } fclose(Potential_file); /* Since the Schrodinger equation is linear, the amplitude of the wavefunction will be fixed by normalization. At left we set it small but nonzero. */ const static double Psi_left = 1.0e-3; // left boundary condition const static double Psi_right = 0.0; // right boundary condition double *Psi, *EigenEnergies;// Arrays to hold the results Psi = new double[N_of_Divisions+1]; //N_of_Points = N_of_Divisions+1 EigenEnergies = new double[N_max+1]; Psi[0] = Psi_left; Psi[1] = Psi_left + 1.0e-3; // Add arbitrary small value int N_quantum;//N_quantum is Energy Quantum Number int Nodes_plus; // Number of nodes (+1) in wavefunction double K_square;// Square of wave vector // Initial Eigen-energy search limits double E_lowerLimit = 0.0;// Eigen-energy must be positive double E_upperLimit = 10.0; int End_sign = -1; bool Limits_are_defined = false; double Normalization_coefficient; double E_trial; // MAIN LOOP begins:----------------------------------- for(N_quantum=1; N_quantum <= N_max; ++N_quantum) { // Find the eigen-values for energy. See theorems (4.1) and (4.2). Limits_are_defined = false; while (Limits_are_defined == false) { /* First, determine an upper limit for energy, so that the wavefunction Psi[i] has one node more than physically needed. */ Nodes_plus = 0; E_trial = E_upperLimit; for (int i=2; i <= N_of_Divisions; ++i) { K_square = 2.0*(E_trial - E_pot[i]); // Now use the NUMEROV-equation (4.197) to calculate wavefunction Psi[i] = 2.0*Psi[i-1]*(1.0 - (5.0*h_0*h_0*K_square / 12.0)) /(1.0 + (h_0*h_0*K_square/12.0))-Psi[i-2]; if (Psi[i]*Psi[i-1] < 0) ++Nodes_plus; } /* If one runs into the following condition, the modification of the upper limit was too aggressive. */ if (E_upperLimit < E_lowerLimit) E_upperLimit = MAX(2*E_upperLimit, -2*E_upperLimit); if (Nodes_plus > N_quantum) E_upperLimit *= 0.7; else if (Nodes_plus < N_quantum) E_upperLimit *= 2.0; else Limits_are_defined = true; // At least one node should appear. } // End of the loop: while (Limits_are_defined == false) // Refine the energy by satisfying the right boundary condition. End_sign = -End_sign; while ((E_upperLimit - E_lowerLimit) > Epsilon) { E_trial = (E_upperLimit + E_lowerLimit) / 2.0; for (int i=2; i <= N_of_Divisions; ++i) { // Again eq.(4.197) of the Numerov-algorithm: K_square = 2.0*(E_trial - E_pot[i]); Psi[i] = 2.0*Psi[i-1] * (1.0 - (5.0*h_0*h_0*K_square / 12.0)) /(1.0 + (h_0*h_0*K_square/12.0))-Psi[i-2]; } if (End_sign*Psi[N_of_Divisions] > Psi_right) E_lowerLimit = E_trial; else E_upperLimit = E_trial; } // End of loop: while ((E_upperLimit - E_lowerLimit) > Epsilon) //APPENDIX C. C++ CODE FOR SOLVING THE SCHRÖDINGER EQUATION // Initialization for the next iteration in main loop E_trial = (E_upperLimit+E_lowerLimit)/2; EigenEnergies[N_quantum] = E_trial; E_upperLimit = E_trial; E_lowerLimit = E_trial; // Now find the normalization coefficient double Integral = 0.0; for (int i=1; i <= N_of_Divisions; ++i) { // Simple integration Integral += 0.5*h_0*(Psi[i-1]*Psi[i-1]+Psi[i]*Psi[i]); } Normalization_coefficient = sqrt(1.0/Integral); // Output of normalized dimensionless wave function for (int i=0; i <=N_of_Divisions; ++i) { fprintf(Wavefunction_file, "%16.12e \t\t %16.12e\n", i*h_0 - 0.5*xRange, Normalization_coefficient*Psi[i]); } fprintf(Wavefunction_file,"\n"); } // End of MAIN LOOP. -------------------------------- fclose(Wavefunction_file); /*Finally convert dimensionless units in real units. Note that energy does not depend explicitly on the particle’s mass anymore: hbar = 1.05457e-34;// Planck constant/2pi omega = 5.34e21; // Frequency in 1/s MeV = 1.602176487e-13; // in J The correct normalization would be hbar*omega/MeV = 3.5148461144, but we use the approximation 3.5 for energy-scale as in chap. 4.9 */ const static double Energyscale = 3.5;// in MeV // Output with rescaled dimensions; assign Energy_file printf("Quantum Harmonic Oscillator, program osci.cpp\n"); printf("Energies in MeV:\n"); printf("n \t\t E_n\n"); for (N_quantum=1; N_quantum <= N_max; ++N_quantum) { fprintf(Energy_file,"%d \t\t %16.12e\n", N_quantum-1, Energyscale*EigenEnergies[N_quantum]); printf("%d \t\t %16.12e\n", N_quantum-1, Energyscale*EigenEnergies[N_quantum]); } fprintf(Energy_file,"\n"); fclose(Energy_file); printf("Wave-Functions in File: psi_osc.dat \n"); printf("\n"); return 0;}

The subnet calculator helps you design networks by taking an IP address and network mask and returning information such as network, broadcast address, and host range.

One of a set of tools we're offering as a way of saying thank you for being a part of the community.

thank you but I'm still getting the following error while compiling with cygwin

DELL@DELL2 ~
$ ./osci.cpp
./osci.cpp: line 1: /bin: Is a directory
./osci.cpp: line 2: a: command not found
./osci.cpp: line 3: To: command not found
./osci.cpp: line 4: To: command not found
./osci.cpp: line 5: Plot: command not found
./osci.cpp: line 6: /GNUPLOT: No such file or directory
./osci.cpp: line 11: $'\r': command not found
./osci.cpp: line 12: syntax error near unexpected token `('
'/osci.cpp: line 12: `int main(int argc, char*argv[])

Here is one way to compile and build an executable named a.exe in Cygwin:

$ g++ osci.cpp

To run the program:
$ ./a.exe
Quantum Harmonic Oscillator, program osci.cpp
Energies in MeV:
n E_n
0 1.749999999795e+00
1 5.249999998112e+00
2 8.749999992829e+00
3 1.224999998232e+01
4 1.574999996759e+01
Wave-Functions in File: psi_osc.dat

/* osci.cpp: Solution of the one-dimensional Schrodinger equation for
a particle in a harmonic potential, using the shooting method.
To compile and link with gnu compiler, type: g++ -o osci osci.cpp
To run the current C++ program, simply type: osci
Plot by gnuplot: /GNUPLOT> set terminal windows
/GNUPLOT> plot "psi-osc.dat" with lines */

the above lines are between /* ... */ what means that the lines are comments. so they should be ignored both by the precompiler and the compiler.

i don't know any c or c++ compiler which would try to compile these comment lines.

did you try the build command suggested by phoffric?

note, at Windows you may consider to using visual studio rather than cygwin. especially if the code is not from unix and/or doesn't use unix specials like fork.

you can download visual studio 2015 community for free. then create a new empty c++ win32 console project (with precompiled headers disabled). add the osci.cpp (i posted) to the project and you are done.

You got those errors because you are trying to run a file that is source code.
C++ is not an interpretive language.

You have to build it as I showed you. After building it using g++, take a look in the same directory and see whether you have the executable named a.exe and if you do, then you can run it by typing
./a.exe at the dollar sign prompt. Note that there is a dot at the beginning of the command followed by a forward slash followed by the name of the executable, in this case the default name, a.exe.

The code compiles without problems after applying Zoppo's correction.

The Questioner has problems with their developer environment which were addressed by Phoffric and answered in further questions the Author has asked for.

Sara

0

Featured Post

The subnet calculator helps you design networks by taking an IP address and network mask and returning information such as network, broadcast address, and host range.

One of a set of tools we're offering as a way of saying thank you for being a part of the community.

Written by John Humphreys
C++ Threading and the POSIX Library
This article will cover the basic information that you need to know in order to make use of the POSIX threading library available for C and C++ on UNIX and most Linux systems.
[s…

C++ Properties
One feature missing from standard C++ that you will find in many other Object Oriented Programming languages is something called a Property (http://www.experts-exchange.com/Programming/Languages/CPP/A_3912-Object-Properties-in-C.ht…

The goal of the video will be to teach the user the concept of local variables and scope. An example of a locally defined variable will be given as well as an explanation of what scope is in C++. The local variable and concept of scope will be relat…

The viewer will be introduced to the member functions push_back and pop_back of the vector class. The video will teach the difference between the two as well as how to use each one along with its functionality.