Untitled

mail@pastecode.io avatar
unknown
c_cpp
2 years ago
5.4 kB
1
Indexable
Never
#include <stdio.h>
#include <stdlib.h>
#include <math.h>



double** make_matrix(int m, int n) {
    double** a;
    int k;
    a = calloc(m, sizeof(double*));
    for (k = 0; k < m; k += 1)
        a[k] = calloc(n, sizeof(double));
    return a;
}

void freeMatrix(double** a,int m) {
    for(int i = 0; i < m; i +=1){
        free(a[i]);
    }
}
double epsilon = pow(10,-9);
typedef struct simplex_t simplex_t;

struct simplex_t {
    int         m;
    int         n;
    int*        var;
    double**    a;
    double*     b;
    double*     x;
    double*     c;
    double      y;
};

void print_lin_program(int m, int n, double * b, double * c, double ** a) {
    int i,j;
    printf("max z = %2.3lfx%d ", c[0], 0);
    for (i = 1; i < n; i +=1) {
        printf("+%2.3lfx%d ", c[i], i);
    }
    printf("\n");
    for (i = 0; i < m; i +=1) {
        for (j = 0; j < n; j+=1) {
            if(j > 0) {
                printf("+ %2.3lfx%d ", a[i][j], j);
            } else {
                printf("%2.3lfx%d ", a[i][j], j);
            }
        }
        printf(" \u2264 %lf \n", b[i]);
    }
}

int init(simplex_t* s, int m, int n, int* var, double** a, double* b, double* x, double* c, double y) {
    int i, k;
    s->m = m;
    s->n = n;
    s->var = var;
    s->a = a;
    s->b = b;
    s->x = x;
    s->c = c;
    s->y = y;
    if(s->var == NULL) {
        s->var = calloc(m+n+1, sizeof(int));
        for(i = 0; i < m+n; i += 1) {
            s->var[i] = i;
        }
    }
    for(k = 0, i = 1; i < m; i += 1) {
        if(b[i] < b[k]) {
            k = i;
        }
    }
    return k;
}

int select_nonbasic(simplex_t* s) {
    int i;
    for(i = 0; i < s->n; i += 1) {
        if(s->c[i] > epsilon) {
            return i;
        }
    }
    return -1;
}

int initial(simplex_t* s, int m, int n, double** a, double* b, double* c, double* x, double y, int* var) {
    int i, j, k;
    double w;
    k = init(s, m, n, var, a, b, x, c, y);
    if(b[k] >= 0) {
        return 1;
    }
}

void pivot(simplex_t* s, int row, int col) {
    double** a = s->a;
    double* b = s->b;
    double* c = s->c;
    int m = s->m;
    int n = s->n;
    int i, j, t;
    t = s->var[col];
    s->var[col] = s->var[n+row];
    s->var[n+row] = t;
    s->y = s->y + c[col] * b[row] / a[row][col];
    for(i = 0; i < n; i += 1) {
        if(i != col) {
            c[i] = c[i] - c[col] * a[row][i] / a[row][col];
        }
    }
    c[col] = - c[col] / a[row][col];
    for(i = 0; i < m; i += 1) {
        if(i != row) {
            b[i] = b[i] - a[i][col] * b[row] / a[row][col];
        }
    }
    for(i = 0; i < m; i += 1) {
        if(i != row) {
            for(j = 0; j < n; j += 1) {
                if(j != col) {
                    a[i][j] = a[i][j] - a[i][col] * a[row][j] / a[row][col];
                }
            }
        }
    }
    for(i = 0; i < m; i = i + 1) {
        if(i != row) {
            a[i][col] = -a[i][col] / a[row][col];
        }
    }
    for(i = 0; i < n; i += 1) {
        if(i != col) {
            a[row][i] = a[row][i] / a[row][col];
        }
    }
    b[row] = b[row] / a[row][col];
    a[row][col] = 1 / a[row][col];
}

double xsimplex(int m, int n, double** a, double* b, double* c, double* x, double y, int* var, int h) {
    simplex_t* s = calloc(1, sizeof(simplex_t));
    int i, row, col;
    if(!initial(s, m, n, a, b, c, x, y, var)) {
        s->var = NULL;
        return -HUGE_VAL; //return nan 
    }
    while((col = select_nonbasic(s)) >= 0) {
        row = -1;
        for(i = 0; i < m; i += 1) {
            if(a[i][col] > epsilon && (row < 0 || b[i] / a[i][col] < b[row] / a[row][col])) {
                row = i;
            }
        }
        if(row < 0) {
            s->var = NULL;
            return HUGE_VAL; //retuen inf
        }
        pivot(s, row, col);
    }
    if (h == 0) {
        for (i = 0; i < n; i += 1) {
            if(s->var[i] < n) {
                x[s->var[i]] = 0;
            }
        }
        for (i = 0; i < m; i += 1) {
            if(s->var[n+1] < n) {
                x[s->var[n+i]] = s->b[i];
            }
        }
        s->var = NULL;
    } else {
        for (i = 0; i < n; i += 1) {
            x[i] = 0;
        }
        for (i = n; i < n + m; i += 1) {
            x[i] = s->b[i-n];
        }
    }
    return s->y;
}

double simplex(int m, int n, double** a, double* b, double* c, double* x, double y) {
    return xsimplex(m, n, a, b, c , x, y, NULL, 0);
}

int main() {
    int m, n, i, j;
    double * b, * c;

    // & berättar för compiler att m och n ska modifieras (läsas som int?)
    scanf("%d %d", &m, &n);

    double ** a = make_matrix(m,n);
    c = calloc(n, sizeof(double));
    b = calloc(m, sizeof(double));

    for(i = 0; i < n; i+=1) {
        scanf("%lf", &c[i]);
    }
    for (i = 0; i < n; i+=1) {
        for (j = 0; j < m; j+=1) {
            scanf("%lf", &a[i][j]);
        }
    }

    for (i = 0; i < m; i+=1) {
        scanf("%lf", &b[i]);
    }
    print_lin_program(m,n,b,c,a);
    double* x = calloc(n+1, sizeof(double));
    double y = 0;
    printf("%lf \n", simplex(m, n, a, b, c, x, y));
    freeMatrix(a,m);
    free(a);
    free(b);
    free(c);
    return 0;

}