CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > General Forums > Main CFD Forum

Best coding practice for c++ solvers

Register Blogs Community New Posts Updated Threads Search

Like Tree12Likes
  • 2 Post By Joachim
  • 1 Post By harishg
  • 3 Post By mprinkey
  • 2 Post By aerosayan
  • 1 Post By aerosayan
  • 3 Post By sbaffini

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   February 23, 2015, 12:49
Default Best coding practice for c++ solvers
  #1
Senior Member
 
Joachim
Join Date: Mar 2012
Location: Paris, France
Posts: 145
Rep Power: 15
Joachim is on a distinguished road
Hi everyone,

so far, I have been working mostly with Fortran for my CFD codes, and I would like to switch to c++. I already know the language (pointers, classes, etc), but I don't have a lot of experience developing large c++ codes.

Could you please recommend me some textbooks/tutorials/etc that detail good coding practice? I have seen a few books/lectures which focused on robustness, with lots of checks all the time. This looks really good, but isn't that slowing down the code a whole lot? I am basically looking for a book that would explain how to develop an elegant and fast code (speed similar to Fortran if possible).

Thank you very much for your help!

Joachim
rajibroy and aero_head like this.
Joachim is offline   Reply With Quote

Old   February 25, 2015, 05:24
Default
  #2
Senior Member
 
N/A
Join Date: Mar 2009
Posts: 189
Rep Power: 17
harishg is on a distinguished road
Theoretically, you can write a fast C++ code or slow Fortran code. You can also wrap c++ around your Fortran code like overture does.

I do not think that any book is a right place to learn good coding practice. Also, the comments will differ from book to book and person to person.

Is there a specific reason to switch to C++? Programming languages are not hard to pick-up. However, if you do not use them for a long time, you will end up forgetting the syntax and would need a reference. Best to stick to the language that you know and switch only if you need to.
harishg is offline   Reply With Quote

Old   February 25, 2015, 11:04
Default
  #3
Senior Member
 
Joachim
Join Date: Mar 2012
Location: Paris, France
Posts: 145
Rep Power: 15
Joachim is on a distinguished road
Well, actually I am getting closer to my graduation date, and I know that most companies use C++ rather than Fortran. I will stick to Fortran for my thesis, but I would like to be as ready as possible with C++ when I start looking for jobs.
Joachim is offline   Reply With Quote

Old   February 25, 2015, 11:16
Default
  #4
Senior Member
 
N/A
Join Date: Mar 2009
Posts: 189
Rep Power: 17
harishg is on a distinguished road
They look more into your work and less into programming skills. Many of my colleagues got into positions without writing a single c++ code.

Good luck.
rajibroy likes this.
harishg is offline   Reply With Quote

Old   February 26, 2015, 16:11
Default Mit ocw
  #5
New Member
 
Ankit Rohatgi
Join Date: Sep 2010
Posts: 5
Rep Power: 16
ankitrohatgi is on a distinguished road
You can always read these notes before going to any C++ specific interview: http://ocw.mit.edu/courses/electrica...lecture-notes/
ankitrohatgi is offline   Reply With Quote

Old   February 26, 2015, 23:56
Default
  #6
Senior Member
 
Michael Prinkey
Join Date: Mar 2009
Location: Pittsburgh PA
Posts: 363
Rep Power: 25
mprinkey will become famous soon enough
There are plenty of general style guides. Google's is pretty good:

http://google-styleguide.googlecode..../cppguide.html

That may give you some idea of what large codes use. It is a place to start. In reality, if you get a job writing CFD code, your employer will likely hand you a style guide before they give you repository access.

While CFD has to serve similar goals as general software (maintainability, extensibility, etc), performance for CFD is far more critical than it is for most general applications. I've written a lot of code in FORTRAN, C, and C++. There are performance pitfalls in every language. Some have a similar root (poor memory ordering, conditionals inside inner loops, etc.). Some are unique to a given language. C++ is a rich and flexible language--look at the OpenFOAM solver codes if you need any convincing of that. But, at the same time, the abstraction and elegance can make code difficult to understand and performance optimization difficult, both for you and for the compiler.

For performance, I tend to write C code, even in C++. I am wary of making use of more advanced features without careful benchmarking. This is important for scalar code, but it is absolutely critical for threaded and vectorized code. Extracting good performance from current and future processors will require OpenMP and SIMD code, because individual cores are not getting much faster--they are just getting more numerous and sprouting wider vector pipelines. See the Intel MIC processors or Nvidia Kepler as good examples. Whatever language you use and whatever coding style you apply needs to take threading and vector operations into primary consideration. They will drive your data structures and your coding style.

There are lots of "High Performance C++" webpages and books, but, honestly, you have to feel your way. I used to avoid using operator overloading in lieu of C-style macros similar to what FLUENT uses:
Code:
#define V_V(A,op,B) do{A[0] op B[0]; A[1] op B[1]; A[2] op B[2];}while(0)
to do basic vector operations. I did that for a good reason. Over the years, I tried to make use of operator overloading for vector operations only to find that the resulting code was much slower than macro'ed or hand-coded C. But with recent compilers, I've found that I can use inline functions and operator overloading to get similar (and even better) performance with something more elegant.
Code:
class vector
{
    private:
    scalar xComp,yComp,zComp;
    public:
    //......lots of other code ....
    inline vector operator+ (const vector A)    // V+A
    {
        return vector(xComp+A.xComp,yComp+A.yComp,zComp+A.zComp);
    }
    inline vector operator- (const vector A)    // V-A
    {
        return vector(xComp-A.xComp,yComp-A.yComp,zComp-A.zComp);
    }
    inline vector operator+= (const vector A)   // V+=A
    {
        xComp+=A.xComp;
        yComp+=A.yComp;
        zComp+=A.zComp;
        return *this;
    }
    inline vector operator-= (const vector A)   // V-=A
    {
        xComp-=A.xComp;
        yComp-=A.yComp;
        zComp-=A.zComp;
        return *this;
    }
}
This allows me to get the cleaner vector notation without compromising performance--that same code 10 years ago ran much slower. You can also use code like this to build chunks of objects to feed to vector operations (AVX/AVX2). It will be code like this that may allow your code to make full use of future CPU features without littering the entire codebase with #pragmas.

I don't know if this is helpful, except to serve as a warning. It is possible to dream up many elaborate tricks and syntactic constructs, but don't be surprised if they produce elegant, slow code. Whatever style/data structure/class/function choices you make, you should prototype and benchmark them thoroughly on the target platform, especially in the context of new and evolving CPU architectures. You then get to pick the most elegant fast one and use that.
mprinkey is offline   Reply With Quote

Old   November 19, 2020, 04:28
Default
  #7
Senior Member
 
Sayan Bhattacharjee
Join Date: Mar 2020
Posts: 495
Rep Power: 8
aerosayan is on a distinguished road
I was looking through the forums and found this question from 2012.
This question is honestly going to be relevant for a few more decades at least.


So, I will leave my answer here for the sake of science. Maybe my experience will be of help to the future.


After all, people will still be asking "Is Fortran faster than C++?" even 30 years later.


Short Answer :


- Don't use anything for performance critical code if you don't know EXACTLY what assembly code will be generated once your code is compiled and what will be the expected and benchmarked performance.

Medium Answer :

The main reason why Fortran is said to be faster than C++ in scientific applications is because of how Fortran and C++ handles pointer aliasing.


You can learn a little about pointer aliasing and its impact on performance here : https://stackoverflow.com/questions/...y-calculations


Many may say that both Fortran and C++ compile down to machine code. If you can produce the same machine code, then it wouldn't be any problem to use either Fortran/C++.


True.


However, the effort and expertise required to compile Fortran and C++ to same( or at least similar) machine code will NEVER be the same.


C++ is an amazing language that allows you full control of your performance critical code. However, if you mess up, everything will go wrong. Just look at the macro based indexing shown in the previous answer.


Code:
#define V_V(A,op,B) do{A[0] op B[0]; A[1] op B[1]; A[2] op B[2];}while(0)
I mean, what could go wrong with a preprocessor macro that is accessible to every portion of your codebase ?



If you truly want to use C++ for your performance critical code, you need to use the __restrict__ directive that is not present in the C++ standard, but is supported by GCC and Clang. This is recommended in the linked stackoverflow answer.



However, as I said before, the amount of effort and expertise required to compile down to same/similar assembly code as Fortran is quite high.


We'll try to optimize the following C++ code :


Code:
void foo(int *array,int *size,int * value)
{
    for(int i=0;i<*size;++i)
    {
        array[i] = 2 * *value;
    }
}
This is a simple example and you can use the __restrict__ directive on the "value" variable to optimize the code. Cool.


What about the "size" variable though?



If you add the __restrict__ directive to the "size" pointer variable, the whole thing gets messed up and introduces multiple hidden branches and jumps in the assembly code. None of which is good for performance.


Again, I'm not saying that you can't obtain good performance in C++, I'm just saying that it would require some heavy effort and expertise to do so. You, and your team really need to know what they are doing.


Anyway, here's the code generated for the example function when __restrict__ is used.


Code:
// URL      : https://gcc.godbolt.org/z/3osMzc
// COMPILER : x86-64 gcc 9.3
// FLAGS    : -c -std=c++14 -O3


// ORIGINAL
void foo(int *array,int *size,int * value)
{
    for(int i=0;i<*size;++i)
    {
        array[i] = 2 * *value;
    }
}

// DISSASSEMBLY
/*
   foo(int*, int*, int*):
        mov     eax, DWORD PTR [rsi]
        test    eax, eax
        jle     .L1
        xor     eax, eax
.L3:
        mov     ecx, DWORD PTR [rdx]       <- UNOPTIMIZED!
        add     ecx, ecx
        mov     DWORD PTR [rdi+rax*4], ecx
        add     rax, 1
        cmp     DWORD PTR [rsi], eax
        jg      .L3
.L1:
        ret
*/

// OPTIMIZED
void bar(int *array,int *size,int * __restrict__ value)
{
    for(int i=0;i<*size;++i)
    {
        array[i] = 2 * *value;
    }
}

// DISSASSEMBLY
/*
bar(int*, int*, int*):
        mov     eax, DWORD PTR [rsi]
        test    eax, eax
        jle     .L1
        mov     edx, DWORD PTR [rdx]        <- OPTIMIZED OUT OF THE LOOP!
        xor     eax, eax
        add     edx, edx
.L3:
        mov     DWORD PTR [rdi+rax*4], edx
        add     rax, 1
        cmp     DWORD PTR [rsi], eax
        jg      .L3
.L1:
        ret
*/

// HORRIBLE
void jar(int *array,int * __restrict__ size,int * __restrict__ value)
{
    for(int i=0;i<*size;++i)
    {
        array[i] = 2 * *value;
    }
}

// DISSASSEMBLY
/*
jar(int*, int*, int*):
        mov     ecx, DWORD PTR [rsi]
        test    ecx, ecx
        jle     .L10
        mov     esi, DWORD PTR [rdx]
        lea     eax, [rcx-1]
        add     esi, esi
        cmp     eax, 2
        jbe     .L15
        mov     edx, ecx
        movd    xmm1, esi
        mov     rax, rdi
        shr     edx, 2
        pshufd  xmm0, xmm1, 0
        sal     rdx, 4
        add     rdx, rdi
.L13:
        movups  XMMWORD PTR [rax], xmm0
        add     rax, 16
        cmp     rax, rdx
        jne     .L13
        mov     eax, ecx
        and     eax, -4
        test    cl, 3
        je      .L17
.L12:
        movsx   rdx, eax
        mov     DWORD PTR [rdi+rdx*4], esi
        lea     edx, [rax+1]
        cmp     edx, ecx
        jge     .L10
        movsx   rdx, edx
        add     eax, 2
        mov     DWORD PTR [rdi+rdx*4], esi
        cmp     ecx, eax
        jle     .L10
        cdqe
        mov     DWORD PTR [rdi+rax*4], esi
.L10:
        ret
.L17:
        ret
.L15:
        xor     eax, eax
        jmp     .L12
*/

/*
int main()
{
    int n = 1000;
    int val = 123;
    int a[n];

    foo(a, &n, &val);
    bar(a, &n, &val);

    return 0;
}
*/
Here's the godbolt link if you wanna play around with it : https://gcc.godbolt.org/z/3osMzc


If you would like to get the maximum performance out of your code, you'll have to be really careful of the assembly code being generated. Even in Fortran : (something I observed earlier) : Fortran->Assembly : How to remove redundant non-vectorized code?




There are other things in C++ which can be problematic :


- Template metaprogramming is awesome if used in moderation. However, if they're overused, the compilation time and executable size skyrocket. Just try to compile OpenFOAM and see how much time it takes. Fortran compilation times are lightning fast.


- There is absolutely no native C++ suitable alternative to Fortran's array operations. In C++, you could use operator overloading, but they're extremely slow. Even though the previous answer says that the code will be inlined, that is only true for simple equations.



- In order to make writing scientific code easier in C++, scientists have actually developed some cool things like Expression Templates. However they HEAVILY depend on C++ Template Metaprogramming using CRTP(Curriously Recurring Template Pattern). This will kill your compilation speeds.


If you still decide to use C++, use a well built library like Eigen3 which uses Expression Templates, but is built to perform well. It just took them huge amount of expertise and effort.




Closing Remarks :


Try giving Fortran a chance.


Many codebases written in Fortran is not extensible, and Fortran gets the bad reputation because developers expect Fortran to do everything.


In my own personal opinion, Fortran codes should be written with the mindset that the code is there just to crunch some numbers as fast as possible and print the results.


Write your code in such a way that it expects the input to be well formatted. Don't just read in the mesh generated by your CAD software. First pass it through a Python/C++ script that renumbers your cells through a Space Filling Curve (like Hilbert Curve/ Z Curve), calculate all of the NodeToCell connectivity information, calculate all of the EdgeToCell connectivity information, calculate all of the NeighbourCell information, order all of your cell nodes in clockwise/counter-clockwise direction, write the renumbered mesh in a simple text file that can be read EASILY by your fortran code.


How I optimized my mesh : How to group cells to distribute onto parallel processors/threads?



Then simply write your Fortran code to do what it does best. Crunch Numbers 4 Life.



Best of luck. Happy coding
ssh123 and billykon like this.

Last edited by aerosayan; November 19, 2020 at 04:30. Reason: Fixed date
aerosayan is offline   Reply With Quote

Old   November 19, 2020, 09:36
Default
  #8
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 2,192
Blog Entries: 29
Rep Power: 39
sbaffini will become famous soon enoughsbaffini will become famous soon enough
Send a message via Skype™ to sbaffini
I politely disagree on the level of confinement you suggest for Fortran. There is no actual reason to avoid reading files or any mesh related stuff. Actually, with modern Fortran, I think there is simply nothing in the solver that can't be done EASILY.

GUIs, hardware and graphics, web, interaction with the OS, libraries in general is where I would suggest switching to another language.

There are certainly some idiosyncracies when using it at large but, do we really want to talk about the ones in C++?
sbaffini is offline   Reply With Quote

Old   November 19, 2020, 11:43
Default
  #9
Senior Member
 
Sayan Bhattacharjee
Join Date: Mar 2020
Posts: 495
Rep Power: 8
aerosayan is on a distinguished road
Quote:
Originally Posted by sbaffini View Post
I politely disagree on the level of confinement you suggest for Fortran. There is no actual reason to avoid reading files or any mesh related stuff. Actually, with modern Fortran, I think there is simply nothing in the solver that can't be done EASILY.

GUIs, hardware and graphics, web, interaction with the OS, libraries in general is where I would suggest switching to another language.

There are certainly some idiosyncracies when using it at large but, do we really want to talk about the ones in C++?

While i agree many things can be done with modern fortran, I still think confining the purpose of the fortran code to only be a number cruncher gives some significant advantages in reducing development time and improving code maintainability.


The main reason I'm biased towards using fortran only for number crunching and C++ for preprocessing and other complex stuff, is due to the benefits of the C++ STL library.


Some of the operations I mentioned eariler require STL containers like std::set, std::map and algorithms such as std::sort. Since fortran doesn't have such a standard library and I don't want to download a third party library that may become outdated one day, it makes it easier to use C++ for those tasks.


During the mesh renumbering task, STL library was instrumental in reducing the development time. The following code shows the step when the boundary cells were pushed to the front of the cell list and the fluid cells were pushed to the end--all using std::sort.



Code:
    // Reorder the cell indices until the marker cells are at the front and
    // the interior cells are at the end of the list. However it should be
    // noted that in both of these sectors, the morton ordering is followed
    // to improve cache performance.
    //
    // Two sorting operations are done to reorder the cells in the desired
    // form.
    //
    for(s32 i=0; i<nmkcells; ++i)
    {
        // Mark all marker cells with -100 to allow sorting
        org[mkcells_found[i]][0] = -100;
    }


    // Organize all of the marker cells at the front
    std::sort(&org[0], &org[0]+g.ncells,
        [](const auto& lhs, const auto& rhs) -> bool
        {
            return lhs.data[0] < rhs.data[0];
        }
        );

    // Organize all of the inner cells at the end
    std::sort(&org[0]+nmkcells, &org[0]+g.ncells,
        [](const auto& lhs, const auto& rhs) -> bool
        {
            return lhs.data[1] < rhs.data[1];
        }
        );

The Node2Cell connectivity and Cell2Cell connectivity calculation required heavy use of std::set and std::map.


If I rewrote all of that in fortran, it would take me more than 3 weeks. In C++ I was able to do it in 3 days--and it is more maintainable and robust.


Now I can focus more time on writing, testing and optimizing my core solver in Fortran.
ssh123 likes this.
aerosayan is offline   Reply With Quote

Old   November 19, 2020, 12:05
Default
  #10
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 2,192
Blog Entries: 29
Rep Power: 39
sbaffini will become famous soon enoughsbaffini will become famous soon enough
Send a message via Skype™ to sbaffini
I see what you mean, but I have few objections:

1) You are actually dealing with integers and that's all you are going to actually deal with, in this sort of stuff. I developed myself, once, my routines for lists, sort, search, etc. and that's it. Fortran doesn't have metaprogramming as C++, but if you're not going to make a library, the effort is very limited and doable. Also, it has advantages (see below).

2) When I developed my routines, as stated above, I could actually tailor them to my real needs (which, in my opinion, is kind of the real thing in Fortran not actually having a STL). For example, when processing grids, I do it in batches of cells instead of one cell at the time. In order to do this, your search and sort routines have lot of peculiarities (in example, to decide if a batch is worth looking into it), in my case even repeated indices. You don't get that with someone else library or a STL.

3) You then end up with both C++ and Fortran in a code and everyone working on it will basically need to know both. Truth is, as you wrote, that you will then need good programmers for the C++ part, while an expert on a certain math/physics topic (say, combustion in CFD) is very likely to not be an expert programmer at all. With Fortran, everyone can just jump in, no matter their programming experience (of course they won't be your code architects in a month)
ssh123, aerosayan and aero_head like this.
sbaffini is offline   Reply With Quote

Reply

Tags
c++, cfd, solver


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Solid Mechanics Solvers added to OpenFOAM Extend bigphil OpenFOAM Announcements from Other Sources 26 October 12, 2017 05:01
FSI solvers - best solvers steven123 OpenFOAM Running, Solving & CFD 0 July 8, 2014 11:26
Possible turbulence modelling bug in SRF solvers otm OpenFOAM Running, Solving & CFD 3 May 29, 2012 05:03
Reaction solvers megacrout OpenFOAM 5 July 8, 2011 11:54
User coding on NT w2k Philip Jones Siemens 0 November 26, 2001 04:43


All times are GMT -4. The time now is 07:23.