This content illustrates how Python can be used in conjunction with C libraries to perform complex scientific computations. It shows the process of defining a C structure in Python using the ctypes library, initializing it with specific parameters, and then passing it to a C library function for computation. The example also demonstrates handling of C arrays and structures from Python, emphasizing interoperability between Python and C for efficient scientific computing.

import ctypes
import numpy as np
from numpy.ctypeslib import ndpointer
import write
 
 
class model_results(ctypes.Structure):
    _fields_=[
        ('time', ctypes.c_double),
        ('csecyr', ctypes.c_double),
        ('modnum', ctypes.c_double),
        ('nshell', ctypes.c_double),
        ('sumass', ctypes.c_double),
        ('teffl', ctypes.c_double),
        ('suluml', ctypes.c_double),
        ('hstot', ctypes.c_double),
        ('ddage', ctypes.c_double),
        ('stllast', ctypes.c_double),
        ('x', ctypes.c_double),
        ('cmixl', ctypes.c_double),
        ('trit', ctypes.c_double),
        ('tril', ctypes.c_double),
        ('ps', ctypes.c_double),
        ('ts', ctypes.c_double),
        ('rs', ctypes.c_double),
        ('z', ctypes.c_double),
        ('smllast', ctypes.c_double),
        ('srllast', ctypes.c_double),
        ('sllast', ctypes.c_double),
        ('spllast', ctypes.c_double),
        ('sdllast', ctypes.c_double),
        ('cgsl', ctypes.c_double),
        ('teffl1', ctypes.c_double),
        ('pn', ctypes.c_double),
        ('tlumx', ctypes.c_double * 7),
        ('cfenv', ctypes.c_double * 3),
        ('elem', ctypes.c_double * 8),
        ('extelem', ctypes.c_double * 4),
        ('sml', ctypes.c_double * 1000),
        ('srl', ctypes.c_double * 1000),
        ('sl', ctypes.c_double * 1000),
        ('spl', ctypes.c_double * 1000),
        ('stl', ctypes.c_double * 1000),
        ('sdl', ctypes.c_double * 1000),
        ('lexcom', ctypes.c_bool),
        ('lrot', ctypes.c_bool),
        ('lexcp', ctypes.c_bool),
        ('ldpree', ctypes.c_bool),
        ('ldiaz', ctypes.c_bool),
        ('lc', ctypes.c_bool),
        ('jcore', ctypes.c_int),
        ('jenv', ctypes.c_int),
        ('mlast', ctypes.c_int),
        ('mlast1', ctypes.c_int)
    ]
 
# defaultArgs = {
#     'sumass': ctypes.c_double(1.0),
#     'teffl1': ctypes.c_double(3.6),
#     'suluml': ctypes.c_double(1.2),
#     'x': ctypes.c_double(0.7),
#     'z': ctypes.c_double(0.02),
#     'pn': ctypes.c_double(1.5),
#     'cmixla': ctypes.c_double(1.26),
#     'ddage': ctypes.c_double(1.00e3),
#     'fmass1': ctypes.c_double(0.00099),
#     'fmass2': ctypes.c_double(1.00000), 'beta': ctypes.c_double(1.00000),
#     'elem':(ctypes.c_double * 8)(*[2.95e-5, 1.73285e-1, 1.8e-3, 5.3152e-2,
#                                    0.0e0,4.82273e-1,0.0e0,0.0e0]),
#     'lexcom': ctypes.c_bool(False)
# }
defaultArgs = {
    'sumass': ctypes.c_double(1.0),
    'teffl1': ctypes.c_double(3.60206),
    'suluml': ctypes.c_double(1.5),
    'x': ctypes.c_double(0.71),
    'z': ctypes.c_double(0.02),
    'pn': ctypes.c_double(1.5),
    'cmixla': ctypes.c_double(1.9126),
    'ddage': ctypes.c_double(1.00e4),
    'fmass1': ctypes.c_double(9.9e-6),
    'fmass2': ctypes.c_double(0.9999),
    'beta': ctypes.c_double(1.00000),
    'elem':(ctypes.c_double * 8)(*[3e-5, 0.0024651398, 2.48862e-05,
        0.00071278576, 0.0, 0.00617626274, 0.0e0,
        0.0e0]),
    'lexcom': ctypes.c_bool(True)
}
 
def run(path, da):
    libnewpoly = ctypes.CDLL(path)
 
    libnewpoly.run.restype = ctypes.c_void_p
 
    test1 = model_results.from_address(libnewpoly.run(
        da['sumass'],
        da['teffl1'],
        da['suluml'],
        da['x'],
        da['z'],
        da['pn'],
        da['cmixla'],
        da['ddage'],
        da['fmass1'],
        da['fmass2'],
        da['beta'],
        da['elem'],
        da['lexcom'])
    )
    write.writeout(
            "test.out",
            test1.time,
            test1.csecyr,
            test1.modnum,
            test1.nshell,
            test1.sumass,
            test1.teffl,
            test1.suluml,
            test1.hstot,
            test1.ddage,
            test1.lexcom,
            test1.lrot,
            test1.lexcp,
            test1.ldpree,
            test1.ldiaz,
            test1.jcore,
            test1.jenv,
            test1.cmixl,
            np.array(test1.tlumx),
            test1.trit,
            test1.tril,
            test1.ps,
            test1.ts,
            test1.rs,
            np.array(test1.cfenv),
            test1.x,
            test1.z,
            np.array(test1.sml),
            np.array(test1.srl),
            np.array(test1.sl),
            np.array(test1.spl),
            np.array(test1.stl),
            np.array(test1.sdl),
            test1.mlast,
            test1.mlast1,
            test1.smllast,
            test1.srllast,
            test1.sllast,
            test1.spllast,
            test1.stllast,
            test1.sdllast,
            test1.lc,
            np.array(test1.elem),
            np.array(test1.extelem),
            test1.cgsl,
            test1.teffl1,
            test1.pn
            )
    return test1
 
if __name__ == "__main__":
    test1 = run('libnewpoly.so', defaultArgs)
    print(np.array(test1.sdl))
 
 
    # print(test1.ckelvin)
    # print(test1.mlast)
    # print(test1.mlast1)
    # print(test1.sml[0:1000])
 
    #
    # libnewpoly.free_model_results(ctypes.byref(test1))
    # del(test1)