#!/bin/bash
# Phase 3: External aerodynamics parametric yaw sweep (solver only)
# Mesh must already be generated. Run blockMesh + snappyHexMesh first.

source /usr/lib/openfoam/openfoam2406/etc/bashrc

BASEDIR="$(cd "$(dirname "$0")" && pwd)"
cd "$BASEDIR"

ANGLES="0 5 10 15 20"

RESULTS="$BASEDIR/results.csv"
echo "yaw_deg,Cd,Cs,Cl,CmYaw" > "$RESULTS"

for ANGLE in $ANGLES; do
    echo ""
    echo "========================================="
    echo "=== Yaw angle: ${ANGLE} degrees ==="
    echo "========================================="

    UX=$(python3 -c "import math; print(f'{10*math.cos(math.radians($ANGLE)):.6f}')")
    UY=$(python3 -c "import math; print(f'{10*math.sin(math.radians($ANGLE)):.6f}')")
    COSX=$(python3 -c "import math; print(f'{math.cos(math.radians($ANGLE)):.6f}')")
    SINX=$(python3 -c "import math; print(f'{math.sin(math.radians($ANGLE)):.6f}')")

    echo "Ux=$UX, Uy=$UY"

    # Clean previous run, keep mesh
    rm -rf [1-9]* processor* postProcessing
    rm -rf 0
    cp -r 0.base 0

    # Write velocity BC with freestream type
    cat > 0/U << UEOF
FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
dimensions      [0 1 -1 0 0 0 0];
internalField   uniform ($UX $UY 0);
boundaryField
{
    inlet
    {
        type            freestream;
        freestreamValue uniform ($UX $UY 0);
    }
    outlet
    {
        type            freestream;
        freestreamValue uniform ($UX $UY 0);
    }
    sides
    {
        type            freestream;
        freestreamValue uniform ($UX $UY 0);
    }
    top
    {
        type            slip;
    }
    bottom
    {
        type            slip;
    }
    droneBody
    {
        type            noSlip;
    }
}
UEOF

    # Write controlDict with correct force directions
    cat > system/controlDict << CEOF
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      controlDict;
}
application     simpleFoam;
startFrom       startTime;
startTime       0;
stopAt          endTime;
endTime         500;
deltaT          1;
writeControl    timeStep;
writeInterval   500;
purgeWrite      1;
writeFormat     ascii;
writePrecision  8;
writeCompression off;
timeFormat      general;
timePrecision   6;
runTimeModifiable true;

functions
{
    forceCoeffs1
    {
        type            forceCoeffs;
        libs            (forces);
        writeControl    timeStep;
        writeInterval   1;
        log             true;
        patches         (droneBody);
        rho             rhoInf;
        rhoInf          1.225;
        liftDir         (0 0 1);
        dragDir         ($COSX $SINX 0);
        CofR            (0 0 0);
        pitchAxis       (0 1 0);
        magUInf         10;
        lRef            0.092;
        Aref            0.006648;
    }

    forces1
    {
        type            forces;
        libs            (forces);
        writeControl    timeStep;
        writeInterval   1;
        log             true;
        patches         (droneBody);
        rho             rhoInf;
        rhoInf          1.225;
        CofR            (0 0 0);
    }
}
CEOF

    echo "=== Running simpleFoam for yaw=$ANGLE ==="
    simpleFoam 2>&1 | tail -20
    echo "--- solver done ---"

    # Extract results
    # Column order: Time Cd Cd(f) Cd(r) Cl Cl(f) Cl(r) CmPitch CmRoll CmYaw Cs Cs(f) Cs(r)
    # Indices:       1    2  3     4     5   6     7      8       9      10   11  12    13
    COEFFS_FILE="postProcessing/forceCoeffs1/0/coefficient.dat"
    if [ ! -f "$COEFFS_FILE" ]; then
        COEFFS_FILE=$(find postProcessing/forceCoeffs1 -name "*.dat" 2>/dev/null | head -1)
    fi

    if [ -n "$COEFFS_FILE" ] && [ -f "$COEFFS_FILE" ]; then
        LAST=$(grep -v "^#" "$COEFFS_FILE" | tail -1)
        NITER=$(grep -v "^#" "$COEFFS_FILE" | wc -l)
        CD=$(echo "$LAST" | awk '{print $2}')
        CS=$(echo "$LAST" | awk '{print $11}')
        CL=$(echo "$LAST" | awk '{print $5}')
        CMY=$(echo "$LAST" | awk '{print $10}')
        echo "Converged in $NITER iterations"
        echo "Results: Cd=$CD, Cs=$CS, Cl=$CL, CmYaw=$CMY"
        echo "$ANGLE,$CD,$CS,$CL,$CMY" >> "$RESULTS"
    else
        echo "WARNING: No coefficient file found for yaw=$ANGLE"
        echo "$ANGLE,NaN,NaN,NaN,NaN" >> "$RESULTS"
    fi

    # Archive
    [ -d postProcessing ] && cp -r postProcessing "postProcessing_yaw${ANGLE}"
done

echo ""
echo "=== Sweep complete ==="
cat "$RESULTS"