Calling PhreeqcRM shared library (.so) from Python

[moein]

I have been trying to call PhreeqcRM functions from Python and I have been successful for most of the functions except one, here you can see how I do it (and Im just putting the part of the code related to that function):

from ctypes import *
import numpy as np
cdll.LoadLibrary('libphreeqcrm-3.5.0.so')
libc = CDLL('libphreeqcrm-3.5.0.so')

ic1 = np.zeros(nxyz*7,dtype=np.int32).reshape(nxyz, 7)
ic2 = np.zeros(nxyz*7,dtype=np.int32).reshape(nxyz, 7)
f1 = np.zeros(nxyz*7,dtype=np.float64).reshape(nxyz, 7)

ic1[:,0] = 1       # Solution 1
ic1[:,1] = -1      # Equilibrium phases none
ic1[:,2] = 1       # Exchange 1
ic1[:,3]= -1      # Surface none
ic1[:,4]= -1      # Gas phase none
ic1[:,5]= -1      # Solid solutions none
ic1[:,6]= -1      # Kinetics none
ic2[:,0]= -1      # Solution none
ic2[:,1]= -1      # Equilibrium phases none
ic2[:,2]= -1      # Exchange none
ic2[:,3]= -1      # Surface none
ic2[:,4]= -1      # Gas phase none
ic2[:,5]= -1      # Solid solutions none
ic2[:,6]= -1      # Kinetics none
f1[:,0]= 1.0      # Mixing fraction ic1 Solution
f1[:,1]= 1.0      # Mixing fraction ic1 Equilibrium phases
f1[:,2]= 1.0      # Mixing fraction ic1 Exchange 1
f1[:,3]= 1.0      # Mixing fraction ic1 Surface
f1[:,4]= 1.0      # Mixing fraction ic1 Gas phase
f1[:,5]= 1.0      # Mixing fraction ic1 Solid solutions
f1[:,6]= 1.0      # Mixing fraction ic1 Kinetics


status = libc.RM_InitialPhreeqc2Module(id,ic1.flatten('F').ctypes,ic2.flatten('F').ctypes,f1.flatten('F').ctypes)


here is what I get:

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
ERROR: Initial condition EQUILIBRIUM_PHASES 1 not found.
Initial condition SURFACE 1 not found.
Initial condition EQUILIBRIUM_PHASES 1 not found.
Initial condition SURFACE 1 not found.
Initial condition EQUILIBRIUM_PHASES 1 not found.
Initial condition SURFACE 1 not found.
Initial condition EQUILIBRIUM_PHASES 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.
Initial condition SURFACE 1 not found.

ERROR: PhreeqcRM failed.
ERROR: Processing initial conditions.
ERROR: InitialPhreeqc2Module: Failure in PhreeqcRM

ERROR: PhreeqcRM failed.
ERROR: PhreeqcRM::InitialPhreeqc2Module
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''

it seems that the ic1 array, when its transported to the function, instead of -1, becomes 1, I was wondering if anyone could help me to resolve the issue. in the following you can find the phreeqc input file:

'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
TITLE Example 11.--Transport and cation exchange.
SOLUTION 0  CaCl2
        units            mmol/kgw
        temp             25.0
        pH               7.0     charge
        pe               12.5    O2(g)   -0.68
        Ca               0.6
        Cl               1.2
SOLUTION 1-40  Initial solution for column
        units            mmol/kgw
        temp             25.0
        pH               7.0     charge
        pe               12.5    O2(g)   -0.68
        Na               1.0
        K                0.2
        N(5)             1.2
END
EXCHANGE 1-40
        -equilibrate 1
        X                0.0011
END
SELECTED_OUTPUT
        -file            advect.sel
        -reset           false
#        -step
        -totals          Na Cl K Ca
USER_PUNCH
  -heading  Temperature Pressure Hyd_K
  10 PUNCH TC, PRESSURE
  20 PUNCH CALLBACK(cell_no, 0, "HYDRAULIC_K")
END
'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''



Thanks - Moein

[dlparkhurst]

I'm guessing it is an addressing issue. Looks like index 0 is converted to 1, and so on.

I don't know Python, but try flattening as "C".