zmeri / PC-SAFT

Functions implementing the PC-SAFT equation of state, including association, electrolyte and dipole terms
GNU General Public License v3.0
49 stars 18 forks source link

portability #83

Closed iurisegtovich closed 3 years ago

iurisegtovich commented 3 years ago

hi, I was wondering which is your primary development platform, linux with g++, or windows with msvc.

I am trying to compile on windows using mingw, after a few adjustments to the src code and a hack in setuptools it worked. however it fails some tests: e.g. acetic acid hres in ">>> test_hres(True)" (Acetic acid, liquid: -36024.93847215227 -38924.64 -7.449526900820988 J/mol) which worked on linux with the same adjustments to source code (Acetic acid, vapor: -15393.874135787593 -15393.63 0.0015859533300081408 J/mol)

i will send you the code adjustments in a pull request and the mingw recipe in another comment, if you are interested

the code adjustements were, i think, language standard issues raising errors in my compiler

the setuptools hack was about using libmingw and leaving libmsvcrt out of the linking step

iurisegtovich commented 3 years ago

how to compile in windows using mingw:

there are many alternative ways to setup a mingw toolchain on a windows 7, 8 or 10 machine

i am going with the way summarized as "get git-bash, get anaconda, run "conda init bash" on the git-bash console, run "conda install -c conda-forge m2w64-toolchain_win-64" and "conda install -c conda-forge make"

with this system i ran the following command to build the pcsaft package

python setup.py build_ext --inplace --global-option --compiler=mingw32

with the following hack to setup.py

diff --git a/setup.py b/setup.py
index 6ac32d4..2d93c91 100644
--- a/setup.py
+++ b/setup.py
@@ -8,6 +8,10 @@ from setuptools import setup, Extension, find_packages
 from Cython.Build import cythonize
 import numpy as np

+import distutils.cygwinccompiler
+distutils.cygwinccompiler.get_msvcr = lambda: []
+#https://stackoverflow.com/questions/11182765/how-can-i-build-my-c-extensions-with-mingw-w64-in-python
+
 ext_modules = [
     Extension("pcsaft",
         sources=["pcsaft.pyx"],
zmeri commented 3 years ago

Yeah, I don't think I've ever compiled this code with mingw. As you probably noticed, I merged your pull request to add those changes you made.

What tests exactly failed when compiling with mingw? Can you give me a full list? Then maybe I can get an idea of where the problem might be.

iurisegtovich commented 3 years ago

the examples with acetic acid, water and methanol failed here.

with the following test script (test_cython_hack.py.txt - modified to run all tests and print errors without raising exceptions ) I got the following output (test_cython_output.txt)

test_cython_output.txt test_cython_hacked.py.txt

zmeri commented 3 years ago

Yes, so it's definitely a problem with the association term because all the associating compounds fail, but none of the systems without associating compounds fail. And, since all the tests fail for the associating compounds, that means the error can't be in the XA derivative calculations (unless there are multiple errors). So it seems that the error is either in these lines, these lines, or the findXA function. Since you already have Mingw on your machine, could you investigate those lines and see where the bug is? It seems like there must be a difference between how Mingw and g++ compile something in those lines.

iurisegtovich commented 3 years ago

thanks for the specifics, I am looking into it, so far I wrote the following test in pure cpp hoping to imitate your acetic acid test in python/cython:

it is quite raw as I am not used to coding in c++

#include <stdio.h>
#include <vector>
#include <string>
#include <cmath>
#include "math.h"
#include "../externals/eigen/Eigen/Dense"
using std::vector;
struct add_args {
    vector<double> m;
    vector<double> s;
    vector<double> e;
    vector<double> k_ij;
    vector<double> e_assoc;
    vector<double> vol_a;
    vector<double> dipm;
    vector<double> dip_num;
    vector<double> z;
    double dielc;
    vector<int> assoc_num;
    vector<int> assoc_scheme;
    vector<double> k_hb;
    vector<double> l_ij;
};
double pcsaft_den_cpp(double t, double p, vector<double> x, int phase, add_args &cppargs);
int main( int argc, const char* argv[] )
{
printf( "\nHello World\n\n" );
double t;
t=305.0;
double p;
p=101325.0;
int phase;
phase=0;
vector<double> x(1) ;
x[0]=1.;
add_args  cppargs = add_args();
cppargs.m.push_back(1.3403);
cppargs.s.push_back(3.8582);
cppargs.e.push_back(211.59);
cppargs.vol_a.push_back(0.075550);
cppargs.e_assoc.push_back(3044.4);
cppargs.assoc_num.push_back(2);
cppargs.assoc_scheme.push_back(0);
cppargs.assoc_scheme.push_back(1);
cppargs.assoc_scheme.push_back(1);
cppargs.assoc_scheme.push_back(0);
double ans;
ans = pcsaft_den_cpp(t, p, x,  phase, cppargs) ;
    printf( "%19.16f", ans );
}

and I added to the first lines of the findXA function:

...
vector<double> XA_find(vector<double> XA_guess, vector<double> delta_ij, double den,
    vector<double> x) {
    /**Iterate over this function in order to solve for XA*/
    int num_sites = XA_guess.size();

   std::cout << "The vector elements are : ";
   for(int i=0; i < XA_guess.size(); i++)
   std::cout << XA_guess[i] << ' ';
   std::cout << '\n';
...

the output differs from linux g++ to mingw

linux:

Hello World

The vector elements are : 0.02 0.02 
The vector elements are : 0.51 0.51 
The vector elements are : 0.755 0.755 
The vector elements are : 0.8775 0.8775
...

mingw:

Hello World

The vector elements are : 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are : 0.02 0.02 
...

next, I will look into those lines that you indicated that do calls to this function and it output

iurisegtovich commented 3 years ago

inputs and outputs of first and second calls:

linux g++

The vector elements are XA_guess: 0.02 0.02 
The vector elements are XA: 1 1 
The vector elements are XA_guess: 0.51 0.51 
The vector elements are XA: 1 1 

mingw

The vector elements are XA_guess: 0.02 0.02 
The vector elements are XA: 1 1 
The vector elements are XA_guess: 0.02 0.02 
The vector elements are XA: 0.771172 0.771172 
iurisegtovich commented 3 years ago

problem at the "dif"

code:

        while ((ctr < 100) && (dif > 1e-15)) {
            ctr += 1;

   std::cout << "The vector elements are @415: ";
   for(int i=0; i < XA_old.size(); i++)
   std::cout << XA_old[i] << ' ';
   std::cout << '\n';

            XA = XA_find(XA_old, delta_ij, den, x_assoc);

   std::cout << "The vector elements are @415 re: ";
   for(int i=0; i < XA_old.size(); i++)
   std::cout << XA_old[i] << ' ';
   std::cout << '\n';

            dif = 0.;
            for (int i = 0; i < num_sites; i++) {
                dif += abs(XA[i] - XA_old[i]);
            }
            for (int i = 0; i < num_sites; i++) {
                XA_old[i] = (XA[i] + XA_old[i]) / 2.0;
            }

   std::cout << "The vector elements are @415 nu: ";
   for(int i=0; i < XA_old.size(); i++)
   std::cout << XA_old[i] << ' ';
      std::cout << "dif= " << dif << ' ';
   std::cout << '\n';

        }

log:


(base) segtovichisv@segtovichisv-dell-luos1604:~/Desktop/github-pcsaft-master_DEBUG/PC-SAFT/tests$ head -n12 log_linux.txt

Hello World

The vector elements are : 305 101325 0 1 1.3403 3.8582 211.59 3044.4 0.07555 0 2 0 1 1 0 
The vector elements are @415: 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are XA: 1 1 
The vector elements are @415 re: 0.02 0.02 
The vector elements are @415 nu: 0.51 0.51 dif= 1.96 
The vector elements are @415: 0.51 0.51 
The vector elements are : 0.51 0.51 
The vector elements are XA: 1 1 
(base) segtovichisv@segtovichisv-dell-luos1604:~/Desktop/github-pcsaft-master_DEBUG/PC-SAFT/tests$ head -n12 log_win.txt

Hello World

The vector elements are : 305 101325 0 1 1.3403 3.8582 211.59 3044.4 0.07555 0 2 0 1 1 0 
The vector elements are @415: 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are XA: 1 1 
The vector elements are @415 re: 0.02 0.02 
The vector elements are @415 nu: 0.51 0.51 dif= 0 
The vector elements are @415: 0.02 0.02 
The vector elements are : 0.02 0.02 
The vector elements are XA: 0.771172 0.771172 
iurisegtovich commented 3 years ago

mingw is interpreting abs function diferrently from linux g++

using the following code: std::cout << "\nDEBUG " << XA[i] <<' '<<XA_old[i]<<' '<<XA[i]-XA_old[i]<<' '<<abs(XA[i]-XA_old[i])<<' '<<std::abs(XA[i]-XA_old[i]); just before the line (around line 415)

...
dif += abs(XA[i] - XA_old[i]);
...

linux g++ output: DEBUG 1 0.02 0.98 0.98 0.98 emphasis: DEBUG 1 0.02 0.98 >>>0.98<<< 0.98

mingw g++ output: DEBUG 1 0.02 0.98 0 0.98 emphasis: DEBUG 1 0.02 0.98 >>>0<<< 0.98

zmeri commented 3 years ago

Just to confirm, did that solve the problem for you when using mingw?

iurisegtovich commented 3 years ago

yes, I did not notice any other problems after adding "std::" to those abs calls in the pull request. I think we can close the issue now. thank you for all the replies :)