void cross_real_array(const real_array_t * x,const real_array_t * y, real_array_t* dest)
{
    /* Assert x and y are vectors */
    /* Assert x and y have size 3 */
    /* Assert dest is a vector */
    /* Assert that dest have size 3*/

    real_set(dest, 0, (real_get(x,1) * real_get(y,2)) - (real_get(x,2) * real_get(y,1)));
    real_set(dest, 1, (real_get(x,2) * real_get(y,0)) - (real_get(x,0) * real_get(y,2)));
    real_set(dest, 2, (real_get(x,0) * real_get(y,1)) - (real_get(x,1) * real_get(y,0)));
}

void cross_real_array(const real_array_t * x,const real_array_t * y, real_array_t* dest)
{
    /* Assert x and y are vectors */
    omc_assert_macro((x->ndims == 1) && (x->dim_size[0] == 3));
    /* Assert y is vector of size 3 */
    omc_assert_macro((y->ndims == 1) && (y->dim_size[0] == 3));
    /* Assert dest is vector of size 3 */
    omc_assert_macro((dest->ndims == 1) && (dest->dim_size[0] == 3));

    real_set(dest, 0, (real_get(*x,1) * real_get(*y,2)) - (real_get(*x,2) * real_get(*y,1)));
    real_set(dest, 1, (real_get(*x,2) * real_get(*y,0)) - (real_get(*x,0) * real_get(*y,2)));
    real_set(dest, 2, (real_get(*x,0) * real_get(*y,1)) - (real_get(*x,1) * real_get(*y,0)));
}

void matrix_real_scalar(modelica_real a, real_array_t* dest)
{
    dest->ndims = 2;
    dest->dim_size[0] = 1;
    dest->dim_size[1] = 1;
    real_set(dest, 0, a);
}

void put_real_matrix_element(modelica_real value, int r, int c, real_array_t* dest)
{
    /* Assert that dest hast correct dimension */
    /* Assert that r and c are valid indices */
    real_set(dest, (r * dest->dim_size[1]) + c, value);
    /* printf("Index %d\n",r*dest->dim_size[1]+c); */
}

void array_real_array(real_array_t* dest,int n,real_array_t* first,...)
{
    int i,j,c,m;
    va_list ap;

    real_array_t **elts = (real_array_t**)malloc(sizeof(real_array_t *) * n);
    assert(elts);
    /* collect all array ptrs to simplify traversal.*/
    va_start(ap,first);
    elts[0] = first;
    for (i = 1; i < n; ++i) {
        elts[i] = va_arg(ap, real_array_t*);
    }
    va_end(ap);

    check_base_array_dim_sizes((const base_array_t **)elts,n);

    for (i = 0, c = 0; i < n; ++i) {
        m = base_array_nr_of_elements(elts[i]);
        for (j = 0; j < m; ++j) {
            real_set(dest, c, real_get(elts[i], j));
            c++;
        }
    }
    free(elts);
}

/* function: transpose_real_array
 *
 * Implementation of transpose(A) for matrix A.
 */
void transpose_real_array(const real_array_t * a, real_array_t* dest)
{
    size_t i;
    size_t j;
    /*  size_t k;*/
    size_t n,m;

    if (a->ndims == 1) {
        copy_real_array_data(a,dest);
        return;
    }

    assert(a->ndims==2 && dest->ndims==2);

    n = a->dim_size[0];
    m = a->dim_size[1];

    assert(dest->dim_size[0] == m && dest->dim_size[1] == n);

    for (i = 0; i < n; ++i) {
        for (j = 0; j < m; ++j) {
            real_set(dest, (j * n) + i, real_get(a, (i * m) + j));
        }
    }
}

void mul_real_matrix_product(const real_array_t * a,const real_array_t * b,real_array_t* dest)
{
    modelica_real tmp;
    size_t i_size;
    size_t j_size;
    size_t k_size;
    size_t i;
    size_t j;
    size_t k;

    /* Assert that dest has correct size */
    i_size = dest->dim_size[0];
    j_size = dest->dim_size[1];
    k_size = a->dim_size[1];

    for (i = 0; i < i_size; ++i) {
        for (j = 0; j < j_size; ++j) {
            tmp = 0;
            for (k = 0; k < k_size; ++k) {
                tmp += real_get(a, (i * k_size) + k)*real_get(b, (k * j_size) + j);
            }
            real_set(dest, (i * j_size) + j, tmp);
        }
    }
}

/* Zero based index */
void simple_indexed_assign_real_array1(const real_array_t * source,
               int i1,
               real_array_t* dest)
{
    /* Assert that source has the correct dimension */
    /* Assert that dest has the correct dimension */
    real_set(dest, i1, real_get(source, i1));
}

void range_real_array(modelica_real start, modelica_real stop, modelica_real inc, real_array_t* dest)
{
    int i;
    /* Assert that dest has correct size */
    for (i = 0; i < dest->dim_size[0]; ++i) {
        real_set(dest, i, start + (i * inc));
    }
}

void skew_real_array(const real_array_t * x,real_array_t* dest)
{
    /* Assert x vector*/
    /* Assert x has size 3*/
    /* Assert dest is 3x3*/
    real_set(dest, 0, 0);
    real_set(dest, 1, -real_get(x, 2));
    real_set(dest, 2, real_get(x, 1));
    real_set(dest, 3, real_get(x, 2));
    real_set(dest, 4, 0);
    real_set(dest, 5, -real_get(x, 0));
    real_set(dest, 6, -real_get(x, 1));
    real_set(dest, 7, real_get(x, 0));
    real_set(dest, 8, 0);
}

void mul_scalar_real_array(modelica_real a,const real_array_t * b,real_array_t* dest)
{
    size_t nr_of_elements;
    size_t i;
    /* Assert that dest has correct size*/
    nr_of_elements = base_array_nr_of_elements(b);
    for (i=0; i < nr_of_elements; ++i) {
        real_set(dest, i, a * real_get(b, i));
    }
}

void mul_real_array(const real_array_t *a,const real_array_t *b,real_array_t* dest)
{
  size_t nr_of_elements;
  size_t i;
  /* Assert that a,b have same sizes? */
  nr_of_elements = base_array_nr_of_elements(*a);
  for(i=0; i < nr_of_elements; ++i) {
    real_set(dest, i, real_get(*a, i) * real_get(*b, i));
  }
}

void simple_indexed_assign_real_array2(const real_array_t * source,
               int i1, int i2,
               real_array_t* dest)
{
    size_t index;
    /* Assert that source has correct dimension */
    /* Assert that dest has correct dimension */
    index = (i1 * source->dim_size[1]) + i2;
    real_set(dest, index, real_get(source, index));
}

void fill_real_array(real_array_t* dest,modelica_real s)
{
    size_t nr_of_elements;
    size_t i;

    nr_of_elements = base_array_nr_of_elements(dest);
    for (i = 0; i < nr_of_elements; ++i) {
        real_set(dest, i, s);
    }
}

void mul_real_array_scalar(const real_array_t * a,modelica_real b,real_array_t* dest)
{
    size_t nr_of_elements;
    size_t i;
    /* Assert that dest has correct size*/
    nr_of_elements = base_array_nr_of_elements(a);
    for (i=0; i < nr_of_elements; ++i) {
        real_set(dest, i, real_get(a, i) * b);
    }
}

void division_real_array_scalar(const real_array_t * a,modelica_real b,real_array_t* dest, const char* division_str)
{
    size_t nr_of_elements;
    size_t i;
    /* Assert that dest has correct size*/
    nr_of_elements = base_array_nr_of_elements(a);
    for (i=0; i < nr_of_elements; ++i) {
        real_set(dest, i, DIVISIONNOTIME(real_get(a, i),b,division_str));
    }
}

void cast_integer_array_to_real(const integer_array_t* a, real_array_t* dest)
{
    int els = base_array_nr_of_elements(a);
    int i;
    clone_base_array_spec(a,dest);
    alloc_real_array_data(dest);
    for (i=0; i<els; i++) {
        real_set(dest, i, (modelica_real)integer_get(a,i));
    }
}

void usub_real_array(real_array_t* a)
{
    size_t nr_of_elements, i;

    nr_of_elements = base_array_nr_of_elements(a);
    for(i = 0; i < nr_of_elements; ++i)
    {
        real_set(a, i, -real_get(a, i));
    }
}

/* function: cat_real_array
 *
 * Concatenates n real arrays along the k:th dimension.
 * k is one based
 */
void cat_real_array(int k, real_array_t* dest, int n,
                    real_array_t* first,...)
{
    va_list ap;
    int i, j, r, c;
    int n_sub = 1, n_super = 1;
    int new_k_dim_size = 0;
    real_array_t **elts = (real_array_t**)malloc(sizeof(real_array_t *) * n);

    assert(elts);
    /* collect all array ptrs to simplify traversal.*/
    va_start(ap, first);
    elts[0] = first;

    for (i = 1; i < n; i++) {
        elts[i] = va_arg(ap,real_array_t*);
    }
    va_end(ap);

    /* check dim sizes of all inputs and dest */
    assert(elts[0]->ndims >= k);
    for (i = 0; i < n; i++) {
        assert(dest->ndims == elts[i]->ndims);
        for (j = 0; j < (k - 1); j++) {
            assert(dest->dim_size[j] == elts[i]->dim_size[j]);
        }
        new_k_dim_size += elts[i]->dim_size[k-1];
        for (j = k; j < elts[0]->ndims; j++) {
            assert(dest->dim_size[j] == elts[i]->dim_size[j]);
        }
    }
    assert(dest->dim_size[k-1] == new_k_dim_size);

    /* calculate size of sub and super structure in 1-dim data representation */
    for (i = 0; i < (k - 1); i++) {
        n_super *= elts[0]->dim_size[i];
    }
    for (i = k; i < elts[0]->ndims; i++) {
        n_sub *= elts[0]->dim_size[i];
    }

    /* concatenation along k-th dimension */
    j = 0;
    for (i = 0; i < n_super; i++) {
        for (c = 0; c < n; c++) {
            int n_sub_k = n_sub * elts[c]->dim_size[k-1];
            for (r = 0; r < n_sub_k; r++) {
                real_set(dest, j,
                            real_get(elts[c], r + (i * n_sub_k)));
                j++;
            }
        }
    }
    free(elts);
}

void div_real_array_scalar(const real_array_t * a,modelica_real b,real_array_t* dest)
{
    size_t nr_of_elements;
    size_t i;
    /* Assert that dest has correct size*/
    /* Do we need to check for b=0? */
    nr_of_elements = base_array_nr_of_elements(a);
    for (i=0; i < nr_of_elements; ++i) {
        real_set(dest, i, real_get(a, i)/b);
    }
}