Himmelblau function with simulated annealing

hb.c

// gcc himmelblau_simulated_annealing.c -lm && ./a.out

#include <math.h>
#include <stdio.h>
#include <stdlib.h> // RAND_MAX

double rand01(void)
{
    return rand() / ((double) RAND_MAX);
}

double hb(double x, double y) {
    return pow(((pow(x,2))+y-11),2) + pow((x+(pow(y,2))-7),2);
}

int main(int argc, char *argv[])
{   
    double x0, y0, k, T0;
    x0 = 2;
    y0 = 1;
    k = 0.1f;
    T0 = 1000.0f;
    
    double x1, y1;
    
    int M, N;
    M = 300;
    N = 15;
    
    double xt, yt;   
    double alpha = 0.85;
    double z = hb(x0, y0);
    double ran_x_1, ran_x_2, ran_y_1, ran_y_2;
    double new_value, current_value;
    double form;
           
    for(int i=0; i < M; i++) {
        for(int j=0; j < N; j++) {
            yt = 0;
            xt = 0;
            
            ran_x_1 = rand01();
            ran_x_2 = rand01();
            ran_y_1 = rand01();
            ran_y_2 = rand01();
            
            if(ran_x_1 >= 0.5) {
                x1 = k * ran_x_2;
            } else {
                x1 = -k * ran_x_2;
            }
            
            if(ran_y_1 >= 0.5) {
                y1 = k * ran_y_2;
            } else {
                y1 = -k * ran_y_2;
            }
                
            xt = x0 + x1;
            yt = y0 + y1;
            new_value = hb(xt, yt);
            current_value = hb(x0, y0);
            form = 1.0 / exp((new_value - current_value) / T0);
            
            if(new_value <= current_value) {
                x0 = xt;
                y0 = yt;
            } else if(rand01() <= form) {
                x0 = xt;
                y0 = yt;
            } else {
                x0 = x0;
                y0 = y0;
            }
        }
        
        z = current_value;
        T0 = alpha * T0;
    }
    
    printf("X is: %f\n", x0);
    printf("Y is: %f\n", y0);
    printf("Minimized value is: %f\n", z);
    return 0;
}