function run_hybrid_worker(batchParams, batchIdx)
    % Set up local cluster for parallel pool
    cluster  = parcluster('local');
    nprocs   = str2double(getenv('SLURM_CPUS_PER_TASK'));
    if isnan(nprocs), nprocs = feature('numcores'); end

    tmpdir = fullfile(getenv('TMPDIR'), sprintf('matlab_job_%d', batchIdx));
    if ~exist(tmpdir, 'dir'); mkdir(tmpdir); end
    cluster.JobStorageLocation = tmpdir;

    pool = parpool(cluster, nprocs);

    nJobs = size(batchParams, 1);
    parfor k = 1:nJobs
        % Unpack parameter tuple
        a_s         = batchParams(k, 1);
        theta_deg   = batchParams(k, 2);
        phi_deg     = batchParams(k, 3);
        N_atoms     = batchParams(k, 4);

        theta_rad   = deg2rad(theta_deg);
        phi_rad     = deg2rad(phi_deg);

        % Create unique save directory
        jobName = sprintf('aS_%03d_theta_%03d_phi_%03d_N_%d', a_s, theta_deg, phi_deg, N_atoms);
        saveDir = fullfile('./Results/Data_3D/GradientDescent', jobName);
        if ~exist(saveDir, 'dir')
            mkdir(saveDir);
        end

        % Options for this run
        OptionsStruct = struct;
        OptionsStruct.NumberOfAtoms           = N_atoms;
        OptionsStruct.DipolarPolarAngle       = theta_rad;
        OptionsStruct.DipolarAzimuthAngle     = phi_rad;
        OptionsStruct.ScatteringLength        = a_s;

        OptionsStruct.TrapFrequencies         = [50, 20, 150];
        OptionsStruct.TrapPotentialType       = 'Harmonic';

        OptionsStruct.NumberOfGridPoints      = [128, 256, 128];
        OptionsStruct.Dimensions              = [30, 50, 30];
        OptionsStruct.UseApproximationForLHY  = true;
        OptionsStruct.IncludeDDICutOff        = true;
        OptionsStruct.CutoffType              = 'Cylindrical';
        OptionsStruct.SimulationMode          = 'EnergyMinimization';
        OptionsStruct.GradientDescentMethod   = 'NonLinearCGD';
        OptionsStruct.MaxIterationsForGD      = 15000;
        OptionsStruct.NoiseScaleFactor        = 0.010;

        OptionsStruct.PlotLive                = false;
        OptionsStruct.JobNumber               = 0;
        OptionsStruct.RunOnGPU                = true;
        OptionsStruct.SaveData                = true;
        OptionsStruct.SaveDirectory           = saveDir;

        options = Helper.convertstruct2cell(OptionsStruct);

        sim = Simulator.DipolarGas(options{:});
        pot = Simulator.Potentials(options{:});
        sim.Potential = pot.trap();

        NumberOfOutputs = 5;
        [~, ~, ~, ~, ~, stats] = Helper.runWithProfiling(@() sim.run(), NumberOfOutputs, saveDir);

        fprintf('Batch %d | Job %d: a_s = %d, theta = %d°, phi = %d°, N = %d | Time = %.2f s\n', ...
            batchIdx, k, a_s, theta_deg, phi_deg, N_atoms, stats.runtime);
    end

    delete(pool);
end