C/C++ Why are my mesh points changing in my 1D linear convection C++ program?

AI Thread Summary
The discussion revolves around issues in a C++ program for simulating 1D linear convection. The main problem identified is that the mesh points are changing unexpectedly due to incorrect initialization of the mesh size variable `dx`, which should be defined as `double dx=2.0/(nx-1);` to avoid integer division. Additionally, there is an out-of-bounds access in the array `u` when using `u[n+1][i]`, as the maximum index exceeds the allocated size. The user is advised to use an integer type for the number of mesh points to prevent confusion and potential errors. Overall, addressing these issues is crucial for the program to function correctly without altering unintended values.
ssatonreb
Messages
5
Reaction score
0
I'm learning c++ and currently trying to write a little programme for 1D linear convection (wave equation).
I have managed set all boundary and initial conditions as well as a mesh.
But I have reached the point where I can't understand why the program does what it does.

In red I highlighted place where I am confused.

If someone could help me please.
In particularly I do not understand why does no matter what I do in the highlighted area of the code my mesh points are changed afterwards.

Code:
//  1D linear convection
#include <iostream>
#include <cmath>
using namespace std;

int main(){
    double nx=10;        // mesh points
    double nt=2;        // time steps 
    double dt=0.01;      // time step size
    double c=1;          // constant (wave propagation velocity)
    double dx=2/(nx-1);  // mesh size    
    int i,n ;
    
    const int NX = nx, NT = nt;     // array size assignment
    double x[NX], u[NT][NX];        // array declaration
    
    memset( u, 0, sizeof u );       // Zeroing out array u[NT][NX]
    memset( x, 0, sizeof x );       // Zeroing out array x[NX]
        
         // Construction of the mesh 
    for (i = 1 ; i < NX ; i++){ 
        x[i]=x[i-1]+dx;
         
    }
        // Initial and Boundary Conditions
    for (i = 0 ; i < NX ; i++){
        if ((x[i]>=0.5) && (x[i]<=1)){
                  u[0][i]=2;
        }
        else
            u[0][i]=1;
    }

// The problem lies here 
// --------------------------------------------------------------------------

    for (n=0; n<NT ; n++){
        for (i=1; i<NX-1; i++) {
            u[n+1][i]=u[n][i]-c*dt/dx*(u[n][i]-u[n][i-1]);
        }
    }
// --------------------------------------------------------------------------       
return 0;}

I really appreciate someones help.
 
Technology news on Phys.org
ssatonreb said:
I'm learning c++ and currently trying to write a little programme for 1D linear convection (wave equation).
I have managed set all boundary and initial conditions as well as a mesh.
But I have reached the point where I can't understand why the program does what it does.

In red I highlighted place where I am confused.

If someone could help me please.
In particularly I do not understand why does no matter what I do in the highlighted area of the code my mesh points are changed afterwards.
Change the line that I marked in red.
ssatonreb said:
Code:
//  1D linear convection
#include <iostream>
#include <cmath>
using namespace std;

int main(){
    double nx=10;        // mesh points
    double nt=2;        // time steps 
    double dt=0.01;      // time step size
    double c=1;          // constant (wave propagation velocity)
    double dx=2/(nx-1);  // mesh size[/color]    
    int i,n ;
    
    const int NX = nx, NT = nt;     // array size assignment
    double x[NX], u[NT][NX];        // array declaration
    
    memset( u, 0, sizeof u );       // Zeroing out array u[NT][NX]
    memset( x, 0, sizeof x );       // Zeroing out array x[NX]
        
         // Construction of the mesh 
    for (i = 1 ; i < NX ; i++){ 
        x[i]=x[i-1]+dx;
         
    }
        // Initial and Boundary Conditions
    for (i = 0 ; i < NX ; i++){
        if ((x[i]>=0.5) && (x[i]<=1)){
                  u[0][i]=2;
        }
        else
            u[0][i]=1;
    }

// The problem lies here 
// --------------------------------------------------------------------------

    for (n=0; n<NT ; n++){
        for (i=1; i<NX-1; i++) {
            u[n+1][i]=u[n][i]-c*dt/dx*(u[n][i]-u[n][i-1]);
        }
    }
// --------------------------------------------------------------------------       
return 0;}

I really appreciate someones help.

The line that I marked in red should be changed to this:
double dx=2.0/(nx-1); // mesh size

What was happening was that dx was being initialized to 0.0, which is not what you wanted. As you had it before, the expression on the right was 2/9, which is 0, due to integer division. That value was being promoted to 0.0 and stored in your dx variable.
 
Mark44 said:
Change the line that I marked in red.


The line that I marked in red should be changed to this:
double dx=2.0/(nx-1); // mesh size

What was happening was that dx was being initialized to 0.0, which is not what you wanted. As you had it before, the expression on the right was 2/9, which is 0, due to integer division. That value was being promoted to 0.0 and stored in your dx variable.

I think this is not a problem:
nx is a double, making (nx-1) a double, making 2/(nx-1) a double.

Another problem however, is that u[n+1] is indexing u out-of-range, which will give unpredictable results.
 
Thank you, Mark44.
I did as you suggested, but it didn't change my original problem, which is related with values of x.
Before I execute lines marked in red values of x are:
0 0.22(2) 0.44(4) 0.66(6) 0.88(8) 1.11(1) 1.33(3) 1.55(5) 1.77(7) 2
and after:
0 0.955 1 1.91203 1.99798 1.08797 1.00202 1 1 2

I have difficulty to understand the reason why does these values change at all because I don't do anything with these values.

Is it related with memory allocation?
 
Sorry - I glanced at nx = 10 and its comment and assumed it was an int. Since you are using it for the number of mesh points, why isn't it an int? When you use a variable to keep track of the number of things, use an integral type, like int or long.

I'm inclined to agree with I like Serena. It seems likely that your statement in red might be tromping over memory that isn't in the array. Your array u runs from u[0][0] through u[1][9].
In your nested loop the largest value of n is NT - 1 == 1 and the largest value of i is NX - 2 == 8.

The statement u[n+1] = ... attempts to store a value in u[2], which is outside this array.
 
I use double for nx because later I use it to calculate dx, which in case of int is just zero.
Thank you very much for your help!
I corrected my mistake and took into account Mark44's first comment I was able also get rid of a line:
const int NX = nx, NT = nt;
Also thank you, I like Serena.
I overlooked this point.
Thank you.
 
Dear Peeps I have posted a few questions about programing on this sectio of the PF forum. I want to ask you veterans how you folks learn program in assembly and about computer architecture for the x86 family. In addition to finish learning C, I am also reading the book From bits to Gates to C and Beyond. In the book, it uses the mini LC3 assembly language. I also have books on assembly programming and computer architecture. The few famous ones i have are Computer Organization and...
What percentage of programmers have learned to touch type? Have you? Do you think it's important, not just for programming, but for more-than-casual computer users generally? ChatGPT didn't have much on it ("Research indicates that less than 20% of people can touch type fluently, with many relying on the hunt-and-peck method for typing ."). 'Hunt-and-peck method' made me smile. It added, "For programmers, touch typing is a valuable skill that can enhance speed, accuracy, and focus. While...
I had a Microsoft Technical interview this past Friday, the question I was asked was this : How do you find the middle value for a dataset that is too big to fit in RAM? I was not able to figure this out during the interview, but I have been look in this all weekend and I read something online that said it can be done at O(N) using something called the counting sort histogram algorithm ( I did not learn that in my advanced data structures and algorithms class). I have watched some youtube...
Back
Top