//#####################################################################
// Copyright 2024, Craig Schroeder, Jing Jin.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
#ifndef __Coordinate_Descent__
#define __Coordinate_Descent__
#include <Core/Nonlinear_Equations/ITERATIVE_SOLVER.h>
#include <functional>

namespace PhysBAM{

template<class T,class TV>
std::tuple<TV,T,int> Coordinate_Descent_Minimum(std::function<T(const TV&)> f, const TV& x0, T init_step, T min_step, int max_iterations)
{
    TV x=x0;
    T y=f(x);
    T h=init_step;
    T inc=2;
    T dec=1/inc;
    int steps_at_size=0;
    int increase_steps=3;
    int it=0;

    for(;it<max_iterations;it++)
    {
        if(it>=max_iterations) break;
        bool progress=false;
        for(int a=0;a<TV::m;a++)
        {
            for(int s=-1;s<2;s+=2)
            {
                while(1)
                {
                    TV w=x.Add_Axis(a,s*h);
                    T z=f(w);
                    if(z>=y) break;
                    progress=true;
                    x=w;
                    y=z;
                    if(++steps_at_size==increase_steps)
                    {
                        h*=inc;
                        steps_at_size=0;
                    }
                }
            }
        }
        if(!progress)
        {
            h*=dec;
            steps_at_size=0;
        }
        if(h<min_step) break;
    }
    return {x,y,it};
}
}
#endif
