Feature/python solver

CMakeLists.txt

@@ -166,9 +166,25 @@ FetchContent_Declare(
     GIT_TAG        "v3.0.1"
 )
 
+# Try to find nlohmann_json from system/conda first, fall back to FetchContent
+find_package(nlohmann_json 3.11 QUIET)
+if(NOT nlohmann_json_FOUND)
+  message(STATUS "nlohmann_json not found, fetching from GitHub")
+  FetchContent_Declare(
+      nlohmann_json
+      GIT_REPOSITORY "https://github.com/nlohmann/json"
+      GIT_TAG        "v3.11.3"
+  )
+else()
+  message(STATUS "Found nlohmann_json: ${nlohmann_json_VERSION}")
+endif()
+
 #set(BUILD_SHARED_LIBS ON)
 set(CMAKE_POSITION_INDEPENDENT_CODE ON)
 FetchContent_MakeAvailable(dylib)
+if(NOT nlohmann_json_FOUND)
+  FetchContent_MakeAvailable(nlohmann_json)
+endif()
 set_target_properties(dylib PROPERTIES UNITY_BUILD OFF)
 
 # Python dependencies (only if building Python bindings)
@@ -321,7 +337,7 @@ endif()
 
 target_link_libraries(simcoon
     PUBLIC ${ARMADILLO_LIBRARIES}
-    PRIVATE dylib
+    PRIVATE dylib nlohmann_json::nlohmann_json
 )
 
 # On Windows with Python bindings, link libsimcoon to carma to use NumPy's memory allocator.
@@ -342,12 +358,9 @@ if(WIN32 AND BUILD_PYTHON_BINDINGS)
   endif()
 endif()
 
-# Compile public executables (with main functions)
-set(EXECUTABLES solver identification L_eff Elastic_props ODF PDF)
-foreach(exe_name ${EXECUTABLES})
-  add_executable(${exe_name} software/${exe_name}.cpp)
-  target_link_libraries(${exe_name} PRIVATE simcoon ${CMAKE_DL_LIBS})
-endforeach()
+# NOTE: Legacy standalone executables (solver, identification, L_eff, Elastic_props, ODF, PDF)
+# have been removed in v2.0. Use the Python API instead for solver workflows.
+# The umat_* source files in software/ are kept for Abaqus/Ansys UMAT builds.
 
 ##Testing
 if(SIMCOON_BUILD_TESTS)

Install.sh

@@ -126,49 +126,12 @@ then
     fi
     Test_OK=$?
 
-    #Create the list of the file to copy after compilation
-    executableToCopy="solver identification L_eff Elastic_props ODF PDF"
-#    objectToCopy="umat_single umat_singleT"
-
     # Copy all important files (+ final message)
+    # NOTE: Legacy executables (solver, identification, L_eff, Elastic_props, ODF, PDF)
+    # have been removed in v2.0. Use the Python API instead.
     if [ $Test_OK -eq 0 ]
     then
         echo "\n---------------------------"
-
-        #Treatement of object files
-#        for object in ${objectToCopy}
-#        do
-#            #Copy of the "object".o from build/CMakeFiles/umat.dir/software to build/bin
-#            if [ -f ${current_dir}/build/CMakeFiles/umat.dir/software/${object}.cpp.o ]
-#            then
-#                cp ${current_dir}/build/CMakeFiles/umat.dir/software/${object}.cpp.o ${current_dir}/build/bin/${object}.o
-#                echo "${blue}${object}.o${reset} copied in ${blue}${current_dir}/build/bin${reset}"
-#            fi
-#        done
-
-        #Treatement of executable files
-        for file in ${executableToCopy}
-        do
-            #if debug exists, copy of the file from build/bin/Debug to build/bin
-            if [ -f ${current_dir}/build/bin/Debug/${file} ]
-            then
-                cp ${current_dir}/build/bin/Debug/${file} ${current_dir}/build/bin
-            fi
-
-            #if Release exists, copy of the file from build/bin/Debug to build/bin
-            if [ -f ${current_dir}/build/bin/Release/${file} ]
-            then
-                cp ${current_dir}/build/bin/Release/${file} ${current_dir}/build/bin
-            fi
-
-            #Copy the file from build/bin to exec
-            cp ${current_dir}/build/bin/${file} ${current_dir}/exec/
-            echo "${blue}${file}${reset} copied in ${blue}${current_dir}/exec${reset}"
-        done
-
-        cp ${current_dir}/build/bin/solver ${current_dir}/examples/elastic-plastic_tension
-        cp ${current_dir}/build/bin/solver ${current_dir}/examples/micromechanics
-        cp ${current_dir}/build/bin/identification ${current_dir}/examples/multi-layer_identification
 
         if [ "${Install_check}" = "OK" ]
         then

README.md

@@ -21,6 +21,44 @@ Simcoon is mainly developed by faculty and researchers from University of Bordea
 
 Simcoon make use and therefore include the FTensor library (http://www.wlandry.net/Projects/FTensor) for convenience. FTensor is a library that handle complex tensor computations. FTensor is released under the GNU General Public License: GPL, version 2. You can get it there (but is is already included in simcoon): (https://bitbucket.org/wlandry/ftensor)
 
+Quick Start (Python)
+--------------------
+
+Simcoon v2.0 introduces a new Python-based solver API that replaces the legacy file-based C++ solver:
+
+```python
+import numpy as np
+from simcoon.solver import Solver, Block, StepMeca
+
+# Material properties: E, nu, alpha
+props = np.array([210000.0, 0.3, 1e-5])
+
+# Define uniaxial tension loading
+step = StepMeca(
+    DEtot_end=np.array([0.01, 0, 0, 0, 0, 0]),  # 1% strain
+    control=['strain', 'stress', 'stress', 'stress', 'stress', 'stress'],
+    Dn_init=50
+)
+
+# Create simulation block
+block = Block(
+    steps=[step],
+    umat_name='ELISO',
+    props=props,
+    nstatev=1
+)
+
+# Run simulation
+solver = Solver(blocks=[block])
+history = solver.solve()
+
+# Extract results
+strain = np.array([h.Etot[0] for h in history])
+stress = np.array([h.sigma[0] for h in history])
+```
+
+For more examples, see `examples/umats/` and the [documentation](https://3mah.github.io/simcoon-docs/).
+
 Documentation
 --------------
 

benchmark/benchmark_comparison.md

@@ -0,0 +1,51 @@
+# Simcoon Solver Benchmark Comparison
+
+## Overview
+
+This benchmark compares the **old C++ file-based solver** (v1.x, master branch) with the
+**new Python solver API** (v2.0, feature/python_solver branch).
+
+## Test Cases
+
+| UMAT | Description |
+|------|-------------|
+| ELISO | Isotropic linear elasticity |
+| EPICP | Plasticity with isotropic hardening |
+| EPKCP | Plasticity with combined isotropic/kinematic hardening |
+| NEOHC | Neo-Hookean hyperelasticity (finite strain) |
+
+## Performance Results
+
+| UMAT | Old Solver (ms) | New Solver (ms) | Speedup | Max Error (%) |
+|------|-----------------|-----------------|---------|---------------|
+| ELISO | 1.75 ± 0.12 | 3.91 ± 0.19 | 0.45x | 0.0000 |
+| EPICP | 1.85 ± 0.05 | 2.69 ± 0.07 | 0.69x | 0.0002 |
+| EPKCP | 2.00 ± 0.10 | 2.81 ± 0.07 | 0.71x | 0.0002 |
+| NEOHC | 3.28 ± 0.14 | 7.54 ± 0.14 | 0.43x | 0.0004 |
+
+## Accuracy Analysis
+
+The maximum relative error between solvers is computed as:
+
+```
+max_error = max(|σ_old - σ_new|) / max(|σ_old|) × 100%
+```
+
+Both solvers should produce nearly identical results (< 0.01% error) for the same
+material model and loading conditions.
+
+## Notes
+
+- **Timing**: Averaged over 10 runs after warmup
+- **Old solver**: Reads/writes configuration files, includes file I/O overhead
+- **New solver**: Pure Python API, no file I/O required
+- **Speedup > 1**: New solver is faster
+- **Speedup < 1**: Old solver is faster
+
+## Conclusion
+
+The old C++ solver is on average **1.8x faster** than the new Python solver.
+
+Maximum relative error across all tests: **0.000381%**
+
+Both solvers produce **numerically equivalent results**.

benchmark/benchmark_new_solver.py

@@ -0,0 +1,130 @@
+"""
+Benchmark for the new Python solver (v2.0).
+Run on feature/python_solver branch.
+
+Saves results to benchmark_results_new.npz
+"""
+
+import numpy as np
+import time
+from pathlib import Path
+
+# Test cases with parameters that work reliably
+BENCHMARK_CASES = {
+    "ELISO": {
+        "props": np.array([210000.0, 0.3, 1e-5]),
+        "nstatev": 1,
+        "strain_max": 0.01,
+        "n_increments": 100,
+        "control_type": "small_strain",
+        "control": ['strain', 'stress', 'stress', 'stress', 'stress', 'stress'],
+        "description": "Isotropic linear elasticity"
+    },
+    "EPICP": {
+        # E, nu, alpha, sigmaY, k, m, kx1, Dx1, kx2, Dx2, kx3, Dx3
+        "props": np.array([70000.0, 0.3, 1e-5, 300.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]),
+        "nstatev": 14,
+        "strain_max": 0.02,
+        "n_increments": 100,
+        "control_type": "small_strain",
+        # Use fully strain-controlled for comparable results
+        "control": ['strain', 'strain', 'strain', 'strain', 'strain', 'strain'],
+        "description": "Plasticity with isotropic hardening"
+    },
+    "EPKCP": {
+        # E, nu, alpha, sigmaY, k, m, kx1, Dx1, kx2, Dx2, kx3, Dx3
+        "props": np.array([70000.0, 0.3, 1e-5, 300.0, 500.0, 1.0, 10000.0, 100.0, 0.0, 0.0, 0.0, 0.0]),
+        "nstatev": 14,
+        "strain_max": 0.03,
+        "n_increments": 100,
+        "control_type": "small_strain",
+        # Use fully strain-controlled for comparable results
+        "control": ['strain', 'strain', 'strain', 'strain', 'strain', 'strain'],
+        "description": "Plasticity with kinematic hardening"
+    },
+    "NEOHC": {
+        # E, nu, alpha (same as ELISO for Neo-Hookean in simcoon)
+        "props": np.array([70000.0, 0.3, 1e-5]),
+        "nstatev": 1,
+        "strain_max": 0.1,
+        "n_increments": 200,
+        "control_type": "logarithmic",
+        # Fully strain-controlled for stability in finite strain
+        "control": ['strain', 'strain', 'strain', 'strain', 'strain', 'strain'],
+        "description": "Neo-Hookean hyperelasticity"
+    },
+}
+
+N_REPEATS = 10
+
+
+def run_benchmark():
+    from simcoon.solver import Solver, Block, StepMeca
+
+    results = {}
+
+    for case_name, cfg in BENCHMARK_CASES.items():
+        print(f"\nBenchmarking {case_name}: {cfg['description']}")
+
+        step = StepMeca(
+            DEtot_end=np.array([cfg["strain_max"], 0, 0, 0, 0, 0]),
+            Dsigma_end=np.zeros(6),
+            control=cfg["control"],
+            Dn_init=cfg["n_increments"],  # Same as old solver
+            Dn_mini=cfg["n_increments"],  # Fixed (no adaptive reduction)
+            Dn_inc=cfg["n_increments"],   # Fixed (no sub-stepping)
+            time=1.0
+        )
+
+        block = Block(
+            steps=[step],
+            umat_name=case_name,
+            props=cfg["props"],
+            nstatev=cfg["nstatev"],
+            control_type=cfg["control_type"]
+        )
+
+        solver = Solver(blocks=[block], max_iter=50, tol=1e-9)
+
+        # Warmup
+        try:
+            history = solver.solve()
+        except Exception as e:
+            print(f"  FAILED: {e}")
+            continue
+
+        # Timing runs
+        times = []
+        for i in range(N_REPEATS):
+            # Re-create solver for fresh state
+            solver = Solver(blocks=[block], max_iter=50, tol=1e-9)
+            start = time.perf_counter()
+            history = solver.solve()
+            elapsed = time.perf_counter() - start
+            times.append(elapsed)
+
+        strain = np.array([h.Etot[0] for h in history])
+        stress = np.array([h.sigma[0] for h in history])
+
+        results[case_name] = {
+            "strain": strain,
+            "stress": stress,
+            "times": np.array(times),
+            "n_points": len(history)
+        }
+
+        print(f"  {len(history)} points, {np.mean(times)*1000:.2f} ± {np.std(times)*1000:.2f} ms")
+
+    # Save results
+    output_file = Path(__file__).parent / "benchmark_results_new.npz"
+    np.savez(output_file, **{f"{k}_{v2}": v1[v2] for k, v1 in results.items() for v2 in v1})
+    print(f"\nResults saved to {output_file}")
+
+    return results
+
+
+if __name__ == "__main__":
+    print("=" * 60)
+    print("Simcoon v2.0 - New Python Solver Benchmark")
+    print("=" * 60)
+    run_benchmark()