decompiled_code:
int __fastcall qdpack_operator_multiply(__int64 a1, __int64 a2, __int64 a3)
{
__int64 v3; // rdi
__int64 v4; // rsi
int result; // eax

if ( a1 && (v3 = *(_QWORD *)(a1 + 16), a2) && v3 && (v4 = *(_QWORD *)(a2 + 16), a3) && v4 && *(_QWORD *)(a3 + 16) )
result = qdpack_matrix_multiply();
else
result = fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_multiply");
return result;
}

canonical_decompiled_code:
int __fastcall func21(__int64 a1, __int64 a2, __int64 a3)
{
__int64 v3; // rdi
__int64 v4; // rsi
int result; // eax

if ( a1 && (v3 = *(_QWORD *)(a1 + 16), a2) && v3 && (v4 = *(_QWORD *)(a2 + 16), a3) && v4 && *(_QWORD *)(a3 + 16) )
result = qdpack_matrix_multiply();
else
result = fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_multiply");
return result;
}

original_code:
void
qdpack_operator_multiply(qdpack_operator_t *A, qdpack_operator_t *B, qdpack_operator_t *C)
{
    if (A && A->data && B && B->data && C && C->data)
    {
        qdpack_matrix_multiply(A->data, B->data, C->data);
    }
    else
    {
        fprintf(stderr, "%s: ERROR: NULL valued operator\n", __PRETTY_FUNCTION__);
    }
}

canonical_original_code:
void qdpack_operator_multiply(struct _qdpack_operator_t *A, struct _qdpack_operator_t *B, struct _qdpack_operator_t *C)
{
    if (A && A->data && B && B->data && C && C->data)
    {
        qdpack_matrix_multiply(A->data, B->data, C->data);
    }
    else
    {
        fprintf(stderr, "%s: ERROR: NULL valued operator\n", __PRETTY_FUNCTION__);
    }
}


repo: jrjohansson/qdpack

decompiled_code:
int __fastcall qdpack_operator_state_multiply(__int64 a1, __int64 a2, __int64 a3)
{
__int64 v3; // r12
double v4; // r14
__int64 v5; // rbp
__int64 v6; // rbx
double v7; // xmm0_8

if ( !a1 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !*(_QWORD *)(a1 + 16) )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !a2 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !*(_QWORD *)(a2 + 16) )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !a3 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
v3 = *(_QWORD *)(a3 + 16);
if ( !v3 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
v4 = gsl_complex_rect(0.0);
v5 = *(_QWORD *)(a2 + 16);
v6 = *(_QWORD *)(a1 + 16);
v7 = gsl_complex_rect(1.0);
return qdpack_matrix_blas_zgemm(111LL, 111LL, v6, v5, v3, v7, 0.0, v4, 0.0);
}

canonical_decompiled_code:
int __fastcall func22(__int64 a1, __int64 a2, __int64 a3)
{
__int64 v3; // r12
double v4; // r14
__int64 v5; // rbp
__int64 v6; // rbx
double v7; // xmm0_8

if ( !a1 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !*(_QWORD *)(a1 + 16) )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !a2 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !*(_QWORD *)(a2 + 16) )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
if ( !a3 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
v3 = *(_QWORD *)(a3 + 16);
if ( !v3 )
return fprintf(stderr, "%s: ERROR: NULL valued operator\n", "qdpack_operator_state_multiply");
v4 = gsl_complex_rect(0.0);
v5 = *(_QWORD *)(a2 + 16);
v6 = *(_QWORD *)(a1 + 16);
v7 = gsl_complex_rect(1.0);
return qdpack_matrix_blas_zgemm(111LL, 111LL, v6, v5, v3, v7, 0.0, v4, 0.0);
}

original_code:
void
qdpack_operator_state_multiply(qdpack_operator_t *A, qdpack_state_t *B, qdpack_state_t *C)
{
    if (A && A->data && B && B->data && C && C->data)
    {
        qdpack_matrix_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), A->data, B->data, (gsl_complex_rect(0.0,0.0)), C->data);
    }
    else
    {
        fprintf(stderr, "%s: ERROR: NULL valued operator\n", __PRETTY_FUNCTION__);
    }
}

canonical_original_code:
void qdpack_operator_state_multiply(struct _qdpack_operator_t *A, struct _qdpack_state *B, struct _qdpack_state *C)
{
    if (A && A->data && B && B->data && C && C->data)
    {
        qdpack_matrix_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), A->data, B->data, (gsl_complex_rect(0.0,0.0)), C->data);
    }
    else
    {
        fprintf(stderr, "%s: ERROR: NULL valued operator\n", __PRETTY_FUNCTION__);
    }
}

variable_types: struct _qdpack_operator_t { int n; struct qdpack_hilbert_space_t { int nsubsys; int *nstates; } *qs; struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *data; } *

variable_types: struct _qdpack_state { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix *data; } *

variable_types:
C
struct _qdpack_state *

return_type: void

user_defined_types: struct _qdpack_operator_t { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix *data; }

user_defined_types: struct _qdpack_matrix { int n; struct gsl_matrix_complex *data; }

user_defined_types: struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct *block; int owner; }


decompiled_code:
double __fastcall qdpack_operator_expectation_value(int *a1, __int64 a2, double a3, double a4)
{
__int64 v4; // r14
__int64 v5; // rbx
double v6; // r12
double v7; // xmm2_8
double v8; // xmm3_8

v4 = qdpack_matrix_alloc(a1[1], *a1);
qdpack_matrix_multiply(*(_QWORD *)(a2 + 16), *((_QWORD *)a1 + 2), v4);
if ( a1[1] <= 0 )
{
v6 = *(double *)&LC0;
}
else
{
v5 = 0LL;
v6 = *(double *)&LC0;
do
{
v7 = qdpack_matrix_get(v4, v5, v5);
v8 = a4;
a4 = *(double *)&LC0;
v6 = gsl_complex_add(v6, 0.0, v7, v8);
++v5;
}
while ( a1[1] > (int)v5 );
}
qdpack_matrix_free(v4);
return v6;
}

canonical_decompiled_code:
double __fastcall func23(int *a1, __int64 a2, double a3, double a4)
{
__int64 v4; // r14
__int64 v5; // rbx
double v6; // r12
double v7; // xmm2_8
double v8; // xmm3_8

v4 = qdpack_matrix_alloc(a1[1], *a1);
qdpack_matrix_multiply(*(_QWORD *)(a2 + 16), *((_QWORD *)a1 + 2), v4);
if ( a1[1] <= 0 )
{
v6 = *(double *)&LC0;
}
else
{
v5 = 0LL;
v6 = *(double *)&LC0;
do
{
v7 = qdpack_matrix_get(v4, v5, v5);
v8 = a4;
a4 = *(double *)&LC0;
v6 = gsl_complex_add(v6, 0.0, v7, v8);
++v5;
}
while ( a1[1] > (int)v5 );
}
qdpack_matrix_free(v4);
return v6;
}

original_code:
gsl_complex
qdpack_operator_expectation_value(qdpack_operator_t *rho, qdpack_operator_t *op)
{
    qdpack_matrix_t *prod;
    gsl_complex tr;
    int i;
    prod = qdpack_matrix_alloc(rho->m, rho->n);
    qdpack_matrix_multiply(op->data, rho->data, prod);
    do {(&tr)->dat[0]=(0.0); (&tr)->dat[1]=(0.0);} while(0);
    for (i = 0; i < rho->m; i++)
    {
        tr = gsl_complex_add(tr, qdpack_matrix_get(prod, i, i));
    }
    qdpack_matrix_free(prod);
    return tr;
}

canonical_original_code:
struct gsl_complex qdpack_operator_expectation_value(struct _qdpack_operator_t *rho, struct _qdpack_operator_t *op)
{
    struct _qdpack_matrix *prod;
    struct gsl_complex tr;
    int i;
    prod = qdpack_matrix_alloc(rho->m, rho->n);
    qdpack_matrix_multiply(op->data, rho->data, prod);
    do {(&tr)->dat[0]=(0.0); (&tr)->dat[1]=(0.0);} while(0);
    for (i = 0; i < rho->m; i++)
    {
        tr = gsl_complex_add(tr, qdpack_matrix_get(prod, i, i));
    }
    qdpack_matrix_free(prod);
    return tr;
}

variable_types:
op
struct _qdpack_operator_t *

variable_types: struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *

variable_types:
tr
struct gsl_complex { double dat[2]; }

variable_types: i int

return_type: struct gsl_complex { double dat[2]; }

user_defined_types: struct qdpack_hilbert_space_t { int nsubsys; int *nstates; }


decompiled_code:
double __fastcall qdpack_operator_trace(__int64 a1, double a2, double a3)
{
__int64 v3; // rbx
double v4; // r12
double v5; // xmm3_8

if ( *(int *)(a1 + 4) <= 0 )
return *(double *)&LC0;
v3 = 0LL;
v4 = *(double *)&LC0;
do
{
qdpack_operator_get(a1, v3, v3);
v5 = a3;
a3 = *(double *)&LC0;
a2 = gsl_complex_add(v4, 0.0, a2, v5);
v4 = a2;
++v3;
}
while ( *(_DWORD *)(a1 + 4) > (int)v3 );
return v4;
}

canonical_decompiled_code:
double __fastcall func24(__int64 a1, double a2, double a3)
{
__int64 v3; // rbx
double v4; // r12
double v5; // xmm3_8

if ( *(int *)(a1 + 4) <= 0 )
return *(double *)&LC0;
v3 = 0LL;
v4 = *(double *)&LC0;
do
{
func12(a1, v3, v3);
v5 = a3;
a3 = *(double *)&LC0;
a2 = gsl_complex_add(v4, 0.0, a2, v5);
v4 = a2;
++v3;
}
while ( *(_DWORD *)(a1 + 4) > (int)v3 );
return v4;
}

original_code:
gsl_complex
qdpack_operator_trace(qdpack_operator_t *rho)
{
    int i;
    gsl_complex tr;
    do {(&tr)->dat[0]=(0.0); (&tr)->dat[1]=(0.0);} while(0);
    for (i = 0; i < rho->m; i++)
    {
        tr = gsl_complex_add(tr, qdpack_operator_get(rho, i, i));
    }
    return tr;
}

canonical_original_code:
struct gsl_complex qdpack_operator_trace(struct _qdpack_operator_t *rho)
{
    int i;
    struct gsl_complex tr;
    do {(&tr)->dat[0]=(0.0); (&tr)->dat[1]=(0.0);} while(0);
    for (i = 0; i < rho->m; i++)
    {
        tr = gsl_complex_add(tr, qdpack_operator_get(rho, i, i));
    }
    return tr;
}

user_defined_types: struct gsl_block_complex_struct { unsigned long size; double *data; }


original_code:
int
qdpack_floquet_modes(qdpack_simulation_t *sim, qdpack_hilbert_space_t *qs)
{
    qdpack_operator_t *U;
    qdpack_matrix_t *Uevec;
    qdpack_matrix_t *Ueval;
    int i, qsn = qdpack_hilbert_space_nstates(qs);
    U = qdpack_operator_alloc(qs);
    qdpack_operator_set_zero(U);
    sim->Ti = 0.0;
    sim->Tf = 2*3.14159265358979323846/sim->h_td_w;
    sim->dT = 2*3.14159265358979323846/sim->h_td_w;
    if (qdpack_propagator(sim, qs, U, ((void *)0)) == -1)
    {
        fprintf(stderr, "Evaluation of the propagator failed.\n");
        return -1;
    }
    Ueval = qdpack_matrix_alloc(qsn, 0);
    Uevec = qdpack_matrix_alloc(qsn, qsn);
    sim->floquet_quasienergies = qdpack_matrix_alloc(qsn, 0);
    sim->floquet_modes = qdpack_operator_alloc(qs);
    qdpack_matrix_eigen_zgeev(U->data, Ueval, Uevec, 3);
    qdpack_matrix_memcpy(sim->floquet_modes->data, Uevec);
    for (i = 0; i < Ueval->m; i++)
    {
        double arg = gsl_complex_arg(qdpack_matrix_get(Ueval, i, 0));
        if (arg > 3.14159265358979323846) arg -= 2*3.14159265358979323846;
        if (arg <= -3.14159265358979323846) arg += 2*3.14159265358979323846;
        qdpack_matrix_set(sim->floquet_quasienergies, i, 0, gsl_complex_rect(-arg/sim->Tf, 0.0));
    }
    qdpack_matrix_free(Uevec);
    qdpack_matrix_free(Ueval);
    qdpack_operator_free(U);
    return 0;
}

canonical_original_code:
int qdpack_floquet_modes(struct _qdpack_simulation *sim, struct qdpack_hilbert_space_t *qs)
{
    struct _qdpack_operator_t *U;
    struct _qdpack_matrix *Uevec;
    struct _qdpack_matrix *Ueval;
    int i; int qsn = qdpack_hilbert_space_nstates(qs);
    U = qdpack_operator_alloc(qs);
    qdpack_operator_set_zero(U);
    sim->Ti = 0.0;
    sim->Tf = 2*3.14159265358979323846/sim->h_td_w;
    sim->dT = 2*3.14159265358979323846/sim->h_td_w;
    if (qdpack_propagator(sim, qs, U, ((void *)0)) == -1)
    {
        fprintf(stderr, "Evaluation of the propagator failed.\n");
        return -1;
    }
    Ueval = qdpack_matrix_alloc(qsn, 0);
    Uevec = qdpack_matrix_alloc(qsn, qsn);
    sim->floquet_quasienergies = qdpack_matrix_alloc(qsn, 0);
    sim->floquet_modes = qdpack_operator_alloc(qs);
    qdpack_matrix_eigen_zgeev(U->data, Ueval, Uevec, 3);
    qdpack_matrix_memcpy(sim->floquet_modes->data, Uevec);
    for (i = 0; i < Ueval->m; i++)
    {
        double arg = gsl_complex_arg(qdpack_matrix_get(Ueval, i, 0));
        if (arg > 3.14159265358979323846) arg -= 2*3.14159265358979323846;
        if (arg <= -3.14159265358979323846) arg += 2*3.14159265358979323846;
        qdpack_matrix_set(sim->floquet_quasienergies, i, 0, gsl_complex_rect(-arg/sim->Tf, 0.0));
    }
    qdpack_matrix_free(Uevec);
    qdpack_matrix_free(Ueval);
    qdpack_operator_free(U);
    return 0;
}

variable_types: struct _qdpack_simulation { struct qdpack_hilbert_space_t { int nsubsys; int *nstates; } *qs; double Ti; double Tf; double dT; char simsig[1024]; double T_w; struct _qdpack_operator_t { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *data; } *H0; struct _qdpack_operator_t *H1; struct _qdpack_operator_t *H_t; struct _qdpack_operator_t *H_a; struct _qdpack_operator_t *H_b; struct _qdpack_operator_t *do_a[20]; struct _qdpack_operator_t *do_ad[20]; struct _qdpack_operator_t *do_aad[20]; struct _qdpack_operator_t *do_ada[20]; struct _qdpack_operator_t *do_a_eb[20]; struct _qdpack_operator_t *do_ad_eb[20]; struct _qdpack_operator_t *do_aad_eb[20]; struct _qdpack_operator_t *do_ada_eb[20]; struct _qdpack_operator_t *ws1; struct _qdpack_operator_t *ws2; struct _qdpack_operator_t *op_tmp; double do_g1[20]; double do_g2[20]; int do_n; void *parameters; struct _qdpack_operator_t *dSSt; double (*ho_w_func)(double, void *); int cont_basis_transform; int option_me_eb; int option_adaptive; double g0; struct _qdpack_operator_t *N_op[10]; struct _qdpack_operator_t *a; struct _qdpack_operator_t *n_op; struct _qdpack_state { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix *data; } *psi_f; double h_td_A; double h_td_w; int (*H0_func)(struct _qdpack_simulation *, struct qdpack_hilbert_space_t *, struct _qdpack_operator_t *); int (*H1_func)(struct _qdpack_simulation *, struct qdpack_hilbert_space_t *, struct _qdpack_operator_t *); int (*Ht_func)(struct _qdpack_simulation *, struct _qdpack_operator_t *, struct _qdpack_operator_t *, struct _qdpack_operator_t *, double); struct _qdpack_operator_t *evec; struct _qdpack_matrix *eval; struct _qdpack_operator_t *floquet_states_t; struct _qdpack_operator_t *floquet_modes; struct _qdpack_matrix *floquet_quasienergies; struct _qdpack_operator_t **floquet_X; struct _qdpack_matrix *floquet_Gamma; struct _qdpack_matrix *floquet_ss_prob; int floquet_k_max; } *

variable_types:
qs
struct qdpack_hilbert_space_t *

variable_types: qsn int

variable_types: arg double

return_type: int

user_defined_types: struct _qdpack_simulation { struct qdpack_hilbert_space_t *qs; double Ti; double Tf; double dT; char simsig[1024]; double T_w; struct _qdpack_operator_t *H0; struct _qdpack_operator_t *H1; struct _qdpack_operator_t *H_t; struct _qdpack_operator_t *H_a; struct _qdpack_operator_t *H_b; struct _qdpack_operator_t *do_a[20]; struct _qdpack_operator_t *do_ad[20]; struct _qdpack_operator_t *do_aad[20]; struct _qdpack_operator_t *do_ada[20]; struct _qdpack_operator_t *do_a_eb[20]; struct _qdpack_operator_t *do_ad_eb[20]; struct _qdpack_operator_t *do_aad_eb[20]; struct _qdpack_operator_t *do_ada_eb[20]; struct _qdpack_operator_t *ws1; struct _qdpack_operator_t *ws2; struct _qdpack_operator_t *op_tmp; double do_g1[20]; double do_g2[20]; int do_n; void *parameters; struct _qdpack_operator_t *dSSt; double (*ho_w_func)(double, void *); int cont_basis_transform; int option_me_eb; int option_adaptive; double g0; struct _qdpack_operator_t *N_op[10]; struct _qdpack_operator_t *a; struct _qdpack_operator_t *n_op; struct _qdpack_state *psi_f; double h_td_A; double h_td_w; int (*H0_func)(struct _qdpack_simulation *, struct qdpack_hilbert_space_t *, struct _qdpack_operator_t *); int (*H1_func)(struct _qdpack_simulation *, struct qdpack_hilbert_space_t *, struct _qdpack_operator_t *); int (*Ht_func)(struct _qdpack_simulation *, struct _qdpack_operator_t *, struct _qdpack_operator_t *, struct _qdpack_operator_t *, double); struct _qdpack_operator_t *evec; struct _qdpack_matrix *eval; struct _qdpack_operator_t *floquet_states_t; struct _qdpack_operator_t *floquet_modes; struct _qdpack_matrix *floquet_quasienergies; struct _qdpack_operator_t **floquet_X; struct _qdpack_matrix *floquet_Gamma; struct _qdpack_matrix *floquet_ss_prob; int floquet_k_max; }


decompiled_code:
__int64 __fastcall qdpack_floquet_modes_t(__int64 a1, __int64 a2, double a3)
{
__int64 v3; // r13
__int64 v4; // r15
double v5; // r12
__int64 v6; // r14
double v7; // xmm1_8
double v8; // xmm0_8
__int64 i; // r12
unsigned int v10; // ebx
double v12; // r15
double v13; // r14
double v14; // xmm1_8
double v15; // xmm2_8
double v16; // xmm3_8
__int64 j; // rbx

v3 = qdpack_operator_alloc(a2);
*(_QWORD *)(a1 + 8) = 0LL;
*(double *)(a1 + 16) = a3;
*(double *)(a1 + 24) = a3;
if ( (unsigned int)qdpack_propagator(a1, a2, v3, 0LL) == -1 )
{
v10 = -1;
fwrite("Evaluation of the propagator failed.\n", 1uLL, 0x25uLL, stderr);
}
else
{
v4 = *(_QWORD *)(a1 + 3144);
v5 = gsl_complex_rect(0.0, 0.0);
v6 = *(_QWORD *)(a1 + 3152);
v7 = 0.0;
v8 = gsl_complex_rect(1.0, 0.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v3, v6, v4, v8, 0.0, v5, 0.0);
for ( i = 0LL; *(_DWORD *)(*(_QWORD *)(a1 + 3160) + 4LL) > (int)i; ++i )
{
for ( j = 0LL; **(_DWORD **)(a1 + 3160) > (int)j; ++j )
{
v12 = qdpack_operator_get(*(_QWORD *)(a1 + 3144), i, j);
v13 = v7;
v14 = qdpack_matrix_get(*(_QWORD *)(a1 + 3168), j, 0LL) * a3;
v15 = gsl_complex_polar(1.0, v14);
v16 = v14;
v7 = v13;
gsl_complex_mul(v12, v13, v15, v16);
qdpack_operator_set(*(_QWORD *)(a1 + 3160), i, j);
}
}
qdpack_operator_free(v3);
v10 = 0;
}
return v10;
}

canonical_decompiled_code:
__int64 __fastcall func1(__int64 a1, __int64 a2, double a3)
{
__int64 v3; // r13
__int64 v4; // r15
double v5; // r12
__int64 v6; // r14
double v7; // xmm1_8
double v8; // xmm0_8
__int64 i; // r12
unsigned int v10; // ebx
double v12; // r15
double v13; // r14
double v14; // xmm1_8
double v15; // xmm2_8
double v16; // xmm3_8
__int64 j; // rbx

v3 = qdpack_operator_alloc(a2);
*(_QWORD *)(a1 + 8) = 0LL;
*(double *)(a1 + 16) = a3;
*(double *)(a1 + 24) = a3;
if ( (unsigned int)qdpack_propagator(a1, a2, v3, 0LL) == -1 )
{
v10 = -1;
fwrite("Evaluation of the propagator failed.\n", 1uLL, 0x25uLL, stderr);
}
else
{
v4 = *(_QWORD *)(a1 + 3144);
v5 = gsl_complex_rect(0.0, 0.0);
v6 = *(_QWORD *)(a1 + 3152);
v7 = 0.0;
v8 = gsl_complex_rect(1.0, 0.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v3, v6, v4, v8, 0.0, v5, 0.0);
for ( i = 0LL; *(_DWORD *)(*(_QWORD *)(a1 + 3160) + 4LL) > (int)i; ++i )
{
for ( j = 0LL; **(_DWORD **)(a1 + 3160) > (int)j; ++j )
{
v12 = qdpack_operator_get(*(_QWORD *)(a1 + 3144), i, j);
v13 = v7;
v14 = qdpack_matrix_get(*(_QWORD *)(a1 + 3168), j, 0LL) * a3;
v15 = gsl_complex_polar(1.0, v14);
v16 = v14;
v7 = v13;
gsl_complex_mul(v12, v13, v15, v16);
qdpack_operator_set(*(_QWORD *)(a1 + 3160), i, j);
}
}
qdpack_operator_free(v3);
v10 = 0;
}
return v10;
}

original_code:
int
qdpack_floquet_modes_t(qdpack_simulation_t *sim,
                       qdpack_hilbert_space_t *qs,
                       double t)
{
    qdpack_operator_t *U;
    int i, j;
    U = qdpack_operator_alloc(qs);
    sim->Ti = 0.0;
       sim->Tf = t;
    sim->dT = t;
    if (qdpack_propagator(sim, qs, U, ((void *)0)) == -1)
    {
        fprintf(stderr, "Evaluation of the propagator failed.\n");
        return -1;
    }
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), U, sim->floquet_modes, (gsl_complex_rect(0.0,0.0)), sim->floquet_states_t);
    for (i = 0; i < sim->floquet_modes_t->m; i++)
    {
        for (j = 0; j < sim->floquet_modes_t->n; j++)
        {
            gsl_complex z1, z2;
            z1 = qdpack_operator_get(sim->floquet_states_t, i, j);
            z2 = gsl_complex_mul(z1, gsl_complex_polar(1.0, ((qdpack_matrix_get(sim->floquet_quasienergies, j, 0)).dat[0]) * t));
            qdpack_operator_set(sim->floquet_modes_t, i, j, z2);
        }
    }
    qdpack_operator_free(U);
    return 0;
}

canonical_original_code:
int qdpack_floquet_modes_t(struct _qdpack_simulation *sim, struct qdpack_hilbert_space_t *qs, double t)
{
    struct _qdpack_operator_t *U;
    int i; int j;
    U = qdpack_operator_alloc(qs);
    sim->Ti = 0.0;
       sim->Tf = t;
    sim->dT = t;
    if (qdpack_propagator(sim, qs, U, ((void *)0)) == -1)
    {
        fprintf(stderr, "Evaluation of the propagator failed.\n");
        return -1;
    }
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), U, sim->floquet_modes, (gsl_complex_rect(0.0,0.0)), sim->floquet_states_t);
    for (i = 0; i < sim->floquet_modes_t->m; i++)
    {
        for (j = 0; j < sim->floquet_modes_t->n; j++)
        {
            struct gsl_complex z1; struct gsl_complex z2;
            z1 = qdpack_operator_get(sim->floquet_states_t, i, j);
            z2 = gsl_complex_mul(z1, gsl_complex_polar(1.0, ((qdpack_matrix_get(sim->floquet_quasienergies, j, 0)).dat[0]) * t));
            qdpack_operator_set(sim->floquet_modes_t, i, j, z2);
        }
    }
    qdpack_operator_free(U);
    return 0;
}

variable_types: t double

variable_types: j int

user_defined_types: struct _qdpack_state { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix *data; }


decompiled_code:
__int64 __fastcall qdpack_floquet_dissipation_rates(__int64 a1, __int64 a2, void (__fastcall *a3)(__int64, __int64, _QWORD, double), void (__fastcall *a4)(__int64, double))
{
void *v4; // rax
__int64 v5; // rbp
_QWORD *v6; // r14
_QWORD *v7; // r14
_QWORD *v8; // r12
double v9; // xmm4_8
int v10; // ebp
__int64 v11; // r15
__int64 v12; // r14
__int64 v13; // r13
double v14; // xmm0_8
double v15; // r13
__int64 v16; // r14
double v17; // xmm0_8
double v18; // xmm1_8
double v19; // xmm0_8
__int64 j; // r13
__int64 k; // rbp
double v22; // r12
double v23; // xmm0_8
__int64 l; // r15
__int64 m; // r14
int v26; // ebp
int v27; // er12
int v28; // er13
__int64 v29; // rbp
double v30; // rbp
double v31; // xmm0_8
double v32; // xmm0_8
double v33; // rbp
double v34; // xmm0_8
double v35; // rbp
double v36; // r12
double v37; // xmm0_8
double v38; // rbp
double v39; // r12
double v40; // xmm1_8
int v41; // er12
__int64 v42; // rbp
__int64 result; // rax
double v44; // [rsp+0h] [rbp-98h]
double v45; // [rsp+0h] [rbp-98h]
double v46; // [rsp+0h] [rbp-98h]
double v47; // [rsp+0h] [rbp-98h]
double v48; // [rsp+0h] [rbp-98h]
double v49; // [rsp+8h] [rbp-90h]
double v50; // [rsp+8h] [rbp-90h]
double v51; // [rsp+10h] [rbp-88h]
double i; // [rsp+18h] [rbp-80h]
int v53; // [rsp+18h] [rbp-80h]
__int64 v54; // [rsp+20h] [rbp-78h]
double v55; // [rsp+28h] [rbp-70h]
__int64 v56; // [rsp+28h] [rbp-70h]
int v57; // [rsp+3Ch] [rbp-5Ch]
double v59; // [rsp+48h] [rbp-50h]
__int64 v60; // [rsp+50h] [rbp-48h]

v57 = qdpack_hilbert_space_nstates(a2);
v54 = qdpack_operator_alloc(a2);
v60 = qdpack_operator_alloc(a2);
*(_DWORD *)(a1 + 3232) = 1;
*(_QWORD *)(a1 + 3184) = malloc(0x18uLL);
*(_QWORD *)(a1 + 3192) = malloc(0x18uLL);
v4 = malloc(0x18uLL);
*(_QWORD *)(a1 + 3176) = v4;
if ( v4 )
{
v5 = 0LL;
do
{
v6 = (_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v5);
*v6 = qdpack_operator_alloc(a2);
qdpack_operator_set_zero(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v5));
v7 = (_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v5);
*v7 = qdpack_operator_alloc(a2);
qdpack_operator_set_zero(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v5));
v8 = (_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v5);
*v8 = qdpack_operator_alloc(a2);
qdpack_operator_set_zero(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v5++));
}
while ( *(_DWORD *)(a1 + 3232) >= (int)v5 - 1 );
*(_QWORD *)(a1 + 3160) = qdpack_operator_alloc(a2);
*(_QWORD *)(a1 + 3144) = qdpack_operator_alloc(a2);
*(_QWORD *)(a1 + 3200) = qdpack_operator_alloc(a2);
if ( !*(_QWORD *)(a1 + 3216) )
*(_QWORD *)(a1 + 3216) = qdpack_matrix_alloc(v57, v57);
a3(a1, a2, *(_QWORD *)(a1 + 3152), 0.0);
v9 = 6.283185307179586 / *(double *)(a1 + 3096);
v49 = v9;
v55 = v9 / 100.0;
v59 = v9 / 100.0 * 0.5 + v9;
for ( i = v9 / 100.0; v59 >= i; i = i + v55 )
{
qdpack_floquet_modes_t(a1, a2, i);
a3(a1, a2, *(_QWORD *)(a1 + 3160), i);
v10 = -*(_DWORD *)(a1 + 3232);
if ( ((*(_DWORD *)(a1 + 3232) & 0x40000000) != 0) == __OFSUB__(*(_DWORD *)(a1 + 3232), v10) )
{
v11 = 0LL;
do
{
v44 = gsl_complex_rect(0.0, 0.0);
v12 = *(_QWORD *)(a1 + 3208);
v13 = *(_QWORD *)(a1 + 3160);
v14 = gsl_complex_rect(1.0, 0.0);
qdpack_operator_blas_zgemm(113LL, 111LL, v13, v12, v54, v14, 0.0, v44, 0.0);
v15 = gsl_complex_rect(0.0, 0.0);
v16 = *(_QWORD *)(a1 + 3160);
v17 = gsl_complex_rect(1.0, 0.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v54, v16, v60, v17, 0.0, v15, 0.0);
v18 = (double)v10 * -6.283185307179586 * i / v49;
v19 = gsl_complex_polar(1.0, v18);
gsl_complex_mul_real(v19, v18, v55 / v49);
qdpack_operator_scale(v60);
qdpack_operator_add(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + v11), v60);
++v10;
v11 += 8LL;
}
while ( *(_DWORD *)(a1 + 3232) >= v10 );
}
}
v53 = -*(_DWORD *)(a1 + 3232);
if ( *(_DWORD *)(a1 + 3232) >= v53 )
{
v56 = 0LL;
do
{
printf("%s: floquet_X[%d]\n", "qdpack_floquet_dissipation_rates", (unsigned int)v56);
qdpack_operator_print_real(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v56));
qdpack_operator_print_imag(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v56));
if ( v57 > 0 )
{
for ( j = 0LL; j != v57; ++j )
{
for ( k = 0LL; k != v57; qdpack_operator_set(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v56), j, k++) )
{
v22 = qdpack_matrix_get(*(_QWORD *)(a1 + 3168), j, 0LL);
v45 = (double)v53 * 6.283185307179586 / v49 + v22 - qdpack_matrix_get(*(_QWORD *)(a1 + 3168), k, 0LL);
qdpack_operator_set(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v56), j, k);
heaviside_function(v45, 0.0);
a4(a1, v45);
v23 = qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v56), j, k);
gsl_complex_abs2(v23);
}
}
}
printf("%s: floquet_delta[%d]\n", "qdpack_floquet_dissipation_rates", (unsigned int)v56);
qdpack_operator_print_real(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v56));
printf("%s: floquet_gamma[%d]\n", "qdpack_floquet_dissipation_rates", (unsigned int)v56);
qdpack_operator_print_real(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v56));
++v53;
++v56;
}
while ( *(_DWORD *)(a1 + 3232) >= v53 );
}
if ( v57 > 0 )
{
for ( l = 0LL; l != v57; ++l )
{
for ( m = 0LL; m != v57; qdpack_operator_set(*(_QWORD *)(a1 + 3200), l, m++) )
{
v26 = *(_DWORD *)(a1 + 3232);
v27 = -v26;
if ( ((v26 & 0x40000000) != 0) == __OFSUB__(v26, -v26) )
{
v46 = 0.0;
do
{
v28 = v26 - v27;
v29 = v27 + v26;
*(_QWORD *)&v50 = COERCE_UNSIGNED_INT64(qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v29), l, m)) & LC14;
v51 = qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v29), l, m);
v30 = qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8LL * v28), m, l);
v46 = (v30 + v51) * boltzmann_factor(v50, *(double *)(a1 + 1056)) + v51 + v46;
++v27;
v26 = *(_DWORD *)(a1 + 3232);
}
while ( v26 >= v27 );
}
}
}
}
printf("%s: floquet_A\n", "qdpack_floquet_dissipation_rates");
qdpack_operator_print_real(*(_QWORD *)(a1 + 3200));
v31 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 0LL, 1LL);
gsl_complex_rect(v31, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v32 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 1LL, 0LL);
gsl_complex_rect(v32, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
v33 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 0LL, 0LL);
v47 = v33 + qdpack_operator_get(*(_QWORD *)(a1 + 3200), 0LL, 1LL);
v48 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 1LL, 0LL) + v47;
v34 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 1LL, 1LL);
gsl_complex_rect((v34 + v48) * 0.5, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 0LL, 0LL);
qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 0LL);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 1LL, 1LL);
if ( !*(_QWORD *)(a1 + 3224) )
*(_QWORD *)(a1 + 3224) = qdpack_matrix_alloc(v57, 0LL);
v35 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
v36 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v37 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
gsl_complex_rect(v35 / (v36 + v37), 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3224), 0LL, 0LL);
v38 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v39 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v40 = v39 + qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
gsl_complex_rect(v38 / v40, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3224), 1LL, 0LL);
v41 = -*(_DWORD *)(a1 + 3232);
if ( !(((*(_DWORD *)(a1 + 3232) & 0x40000000) != 0) ^ __OFSUB__(*(_DWORD *)(a1 + 3232), v41) | (*(_DWORD *)(a1 + 3232) == v41)) )
{
v42 = 0LL;
do
{
qdpack_operator_free(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + v42));
qdpack_operator_free(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + v42));
qdpack_operator_free(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + v42));
++v41;
v42 += 8LL;
}
while ( *(_DWORD *)(a1 + 3232) > v41 );
}
qdpack_operator_free(v54);
qdpack_operator_free(v60);
qdpack_operator_free(*(_QWORD *)(a1 + 3200));
free(*(void **)(a1 + 3176));
free(*(void **)(a1 + 3192));
free(*(void **)(a1 + 3184));
result = 0LL;
}
else
{
fprintf(stderr, "%s: Failed to allocate floquet dissipation matrices.\n", "qdpack_floquet_dissipation_rates");
result = 0xFFFFFFFFLL;
}
return result;
}

canonical_decompiled_code:
__int64 __fastcall func2(__int64 a1, __int64 a2, void (__fastcall *a3)(__int64, __int64, _QWORD, double), void (__fastcall *a4)(__int64, double))
{
void *v4; // rax
__int64 v5; // rbp
_QWORD *v6; // r14
_QWORD *v7; // r14
_QWORD *v8; // r12
double v9; // xmm4_8
int v10; // ebp
__int64 v11; // r15
__int64 v12; // r14
__int64 v13; // r13
double v14; // xmm0_8
double v15; // r13
__int64 v16; // r14
double v17; // xmm0_8
double v18; // xmm1_8
double v19; // xmm0_8
__int64 j; // r13
__int64 k; // rbp
double v22; // r12
double v23; // xmm0_8
__int64 l; // r15
__int64 m; // r14
int v26; // ebp
int v27; // er12
int v28; // er13
__int64 v29; // rbp
double v30; // rbp
double v31; // xmm0_8
double v32; // xmm0_8
double v33; // rbp
double v34; // xmm0_8
double v35; // rbp
double v36; // r12
double v37; // xmm0_8
double v38; // rbp
double v39; // r12
double v40; // xmm1_8
int v41; // er12
__int64 v42; // rbp
__int64 result; // rax
double v44; // [rsp+0h] [rbp-98h]
double v45; // [rsp+0h] [rbp-98h]
double v46; // [rsp+0h] [rbp-98h]
double v47; // [rsp+0h] [rbp-98h]
double v48; // [rsp+0h] [rbp-98h]
double v49; // [rsp+8h] [rbp-90h]
double v50; // [rsp+8h] [rbp-90h]
double v51; // [rsp+10h] [rbp-88h]
double i; // [rsp+18h] [rbp-80h]
int v53; // [rsp+18h] [rbp-80h]
__int64 v54; // [rsp+20h] [rbp-78h]
double v55; // [rsp+28h] [rbp-70h]
__int64 v56; // [rsp+28h] [rbp-70h]
int v57; // [rsp+3Ch] [rbp-5Ch]
double v59; // [rsp+48h] [rbp-50h]
__int64 v60; // [rsp+50h] [rbp-48h]

v57 = qdpack_hilbert_space_nstates(a2);
v54 = qdpack_operator_alloc(a2);
v60 = qdpack_operator_alloc(a2);
*(_DWORD *)(a1 + 3232) = 1;
*(_QWORD *)(a1 + 3184) = malloc(0x18uLL);
*(_QWORD *)(a1 + 3192) = malloc(0x18uLL);
v4 = malloc(0x18uLL);
*(_QWORD *)(a1 + 3176) = v4;
if ( v4 )
{
v5 = 0LL;
do
{
v6 = (_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v5);
*v6 = qdpack_operator_alloc(a2);
qdpack_operator_set_zero(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v5));
v7 = (_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v5);
*v7 = qdpack_operator_alloc(a2);
qdpack_operator_set_zero(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v5));
v8 = (_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v5);
*v8 = qdpack_operator_alloc(a2);
qdpack_operator_set_zero(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v5++));
}
while ( *(_DWORD *)(a1 + 3232) >= (int)v5 - 1 );
*(_QWORD *)(a1 + 3160) = qdpack_operator_alloc(a2);
*(_QWORD *)(a1 + 3144) = qdpack_operator_alloc(a2);
*(_QWORD *)(a1 + 3200) = qdpack_operator_alloc(a2);
if ( !*(_QWORD *)(a1 + 3216) )
*(_QWORD *)(a1 + 3216) = qdpack_matrix_alloc(v57, v57);
a3(a1, a2, *(_QWORD *)(a1 + 3152), 0.0);
v9 = 6.283185307179586 / *(double *)(a1 + 3096);
v49 = v9;
v55 = v9 / 100.0;
v59 = v9 / 100.0 * 0.5 + v9;
for ( i = v9 / 100.0; v59 >= i; i = i + v55 )
{
func1(a1, a2, i);
a3(a1, a2, *(_QWORD *)(a1 + 3160), i);
v10 = -*(_DWORD *)(a1 + 3232);
if ( ((*(_DWORD *)(a1 + 3232) & 0x40000000) != 0) == __OFSUB__(*(_DWORD *)(a1 + 3232), v10) )
{
v11 = 0LL;
do
{
v44 = gsl_complex_rect(0.0, 0.0);
v12 = *(_QWORD *)(a1 + 3208);
v13 = *(_QWORD *)(a1 + 3160);
v14 = gsl_complex_rect(1.0, 0.0);
qdpack_operator_blas_zgemm(113LL, 111LL, v13, v12, v54, v14, 0.0, v44, 0.0);
v15 = gsl_complex_rect(0.0, 0.0);
v16 = *(_QWORD *)(a1 + 3160);
v17 = gsl_complex_rect(1.0, 0.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v54, v16, v60, v17, 0.0, v15, 0.0);
v18 = (double)v10 * -6.283185307179586 * i / v49;
v19 = gsl_complex_polar(1.0, v18);
gsl_complex_mul_real(v19, v18, v55 / v49);
qdpack_operator_scale(v60);
qdpack_operator_add(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + v11), v60);
++v10;
v11 += 8LL;
}
while ( *(_DWORD *)(a1 + 3232) >= v10 );
}
}
v53 = -*(_DWORD *)(a1 + 3232);
if ( *(_DWORD *)(a1 + 3232) >= v53 )
{
v56 = 0LL;
do
{
printf("%s: floquet_X[%d]\n", "qdpack_floquet_dissipation_rates", (unsigned int)v56);
qdpack_operator_print_real(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v56));
qdpack_operator_print_imag(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v56));
if ( v57 > 0 )
{
for ( j = 0LL; j != v57; ++j )
{
for ( k = 0LL; k != v57; qdpack_operator_set(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v56), j, k++) )
{
v22 = qdpack_matrix_get(*(_QWORD *)(a1 + 3168), j, 0LL);
v45 = (double)v53 * 6.283185307179586 / v49 + v22 - qdpack_matrix_get(*(_QWORD *)(a1 + 3168), k, 0LL);
qdpack_operator_set(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v56), j, k);
heaviside_function(v45, 0.0);
a4(a1, v45);
v23 = qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + 8 * v56), j, k);
gsl_complex_abs2(v23);
}
}
}
printf("%s: floquet_delta[%d]\n", "qdpack_floquet_dissipation_rates", (unsigned int)v56);
qdpack_operator_print_real(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v56));
printf("%s: floquet_gamma[%d]\n", "qdpack_floquet_dissipation_rates", (unsigned int)v56);
qdpack_operator_print_real(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v56));
++v53;
++v56;
}
while ( *(_DWORD *)(a1 + 3232) >= v53 );
}
if ( v57 > 0 )
{
for ( l = 0LL; l != v57; ++l )
{
for ( m = 0LL; m != v57; qdpack_operator_set(*(_QWORD *)(a1 + 3200), l, m++) )
{
v26 = *(_DWORD *)(a1 + 3232);
v27 = -v26;
if ( ((v26 & 0x40000000) != 0) == __OFSUB__(v26, -v26) )
{
v46 = 0.0;
do
{
v28 = v26 - v27;
v29 = v27 + v26;
*(_QWORD *)&v50 = COERCE_UNSIGNED_INT64(qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + 8 * v29), l, m)) & LC14;
v51 = qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8 * v29), l, m);
v30 = qdpack_operator_get(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + 8LL * v28), m, l);
v46 = (v30 + v51) * boltzmann_factor(v50, *(double *)(a1 + 1056)) + v51 + v46;
++v27;
v26 = *(_DWORD *)(a1 + 3232);
}
while ( v26 >= v27 );
}
}
}
}
printf("%s: floquet_A\n", "qdpack_floquet_dissipation_rates");
qdpack_operator_print_real(*(_QWORD *)(a1 + 3200));
v31 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 0LL, 1LL);
gsl_complex_rect(v31, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v32 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 1LL, 0LL);
gsl_complex_rect(v32, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
v33 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 0LL, 0LL);
v47 = v33 + qdpack_operator_get(*(_QWORD *)(a1 + 3200), 0LL, 1LL);
v48 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 1LL, 0LL) + v47;
v34 = qdpack_operator_get(*(_QWORD *)(a1 + 3200), 1LL, 1LL);
gsl_complex_rect((v34 + v48) * 0.5, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 0LL, 0LL);
qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 0LL);
qdpack_matrix_set(*(_QWORD *)(a1 + 3216), 1LL, 1LL);
if ( !*(_QWORD *)(a1 + 3224) )
*(_QWORD *)(a1 + 3224) = qdpack_matrix_alloc(v57, 0LL);
v35 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
v36 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v37 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
gsl_complex_rect(v35 / (v36 + v37), 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3224), 0LL, 0LL);
v38 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v39 = qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 0LL, 1LL);
v40 = v39 + qdpack_matrix_get(*(_QWORD *)(a1 + 3216), 1LL, 0LL);
gsl_complex_rect(v38 / v40, 0.0);
qdpack_matrix_set(*(_QWORD *)(a1 + 3224), 1LL, 0LL);
v41 = -*(_DWORD *)(a1 + 3232);
if ( !(((*(_DWORD *)(a1 + 3232) & 0x40000000) != 0) ^ __OFSUB__(*(_DWORD *)(a1 + 3232), v41) | (*(_DWORD *)(a1 + 3232) == v41)) )
{
v42 = 0LL;
do
{
qdpack_operator_free(*(_QWORD *)(*(_QWORD *)(a1 + 3176) + v42));
qdpack_operator_free(*(_QWORD *)(*(_QWORD *)(a1 + 3192) + v42));
qdpack_operator_free(*(_QWORD *)(*(_QWORD *)(a1 + 3184) + v42));
++v41;
v42 += 8LL;
}
while ( *(_DWORD *)(a1 + 3232) > v41 );
}
qdpack_operator_free(v54);
qdpack_operator_free(v60);
qdpack_operator_free(*(_QWORD *)(a1 + 3200));
free(*(void **)(a1 + 3176));
free(*(void **)(a1 + 3192));
free(*(void **)(a1 + 3184));
result = 0LL;
}
else
{
fprintf(stderr, "%s: Failed to allocate floquet dissipation matrices.\n", "qdpack_floquet_dissipation_rates");
result = 0xFFFFFFFFLL;
}
return result;
}

original_code:
int
qdpack_floquet_dissipation_rates(qdpack_simulation_t *sim,
                                       qdpack_hilbert_space_t *qs,
                                       operator_cb_func_t floquet_modes_store_cb,
                                       spectral_density_cb_t s_func_cb)
{
    gsl_complex z;
    double t, delta_t, T;
    int i, i_limit = 100, j, qsn = qdpack_hilbert_space_nstates(qs);
    qdpack_operator_t *w1, *w2;
    int k, k_idx;
    w1 = qdpack_operator_alloc(qs);
    w2 = qdpack_operator_alloc(qs);
    sim->floquet_k_max = 1;
    sim->floquet_gamma = malloc(sizeof(qdpack_operator_t *) * (2 * sim->floquet_k_max + 1));
    sim->floquet_delta = malloc(sizeof(qdpack_operator_t *) * (2 * sim->floquet_k_max + 1));
    sim->floquet_X = malloc(sizeof(qdpack_operator_t *) * (2 * sim->floquet_k_max + 1));
    if (sim->floquet_X == ((void *)0))
    {
        fprintf(stderr, "%s: Failed to allocate floquet dissipation matrices.\n", __PRETTY_FUNCTION__);
        return -1;
    }
    for (k_idx = 0, k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++, k_idx++)
    {
        sim->floquet_X[k_idx] = qdpack_operator_alloc(qs);
        qdpack_operator_set_zero(sim->floquet_X[k_idx]);
        sim->floquet_delta[k_idx] = qdpack_operator_alloc(qs);
        qdpack_operator_set_zero(sim->floquet_delta[k_idx]);
        sim->floquet_gamma[k_idx] = qdpack_operator_alloc(qs);
        qdpack_operator_set_zero(sim->floquet_gamma[k_idx]);
    }
    sim->floquet_modes_t = qdpack_operator_alloc(qs);
    sim->floquet_states_t = qdpack_operator_alloc(qs);
    sim->floquet_A = qdpack_operator_alloc(qs);
    if (sim->floquet_Gamma == ((void *)0))
    {
        sim->floquet_Gamma = qdpack_matrix_alloc(qsn, qsn);
    }
    floquet_modes_store_cb(sim, qs, sim->floquet_modes, 0.0);
    T = 2*3.14159265358979323846/sim->h_td_w;
    delta_t = T / i_limit;
    for (t = delta_t; t <= T+delta_t/2.0; t += delta_t)
    {
        qdpack_floquet_modes_t(sim, qs, t);
        floquet_modes_store_cb(sim, qs, sim->floquet_modes_t, t);
        for (k_idx = 0, k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++, k_idx++)
        {
            qdpack_operator_blas_zgemm(CblasConjTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), sim->floquet_modes_t, sim->floquet_operator, (gsl_complex_rect(0.0,0.0)), w1);
            qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), w1 , sim->floquet_modes_t, (gsl_complex_rect(0.0,0.0)), w2);
            z = gsl_complex_mul_real(gsl_complex_polar(1.0, - 2 * 3.14159265358979323846 * k * t / T), delta_t / T);
            qdpack_operator_scale(w2, z);
            qdpack_operator_add(sim->floquet_X[k_idx], w2);
        }
    }
    for (k_idx = 0, k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++, k_idx++)
    {
        printf("%s: floquet_X[%d]\n", __PRETTY_FUNCTION__, k_idx);
        qdpack_operator_print_real(sim->floquet_X[k_idx]);
        qdpack_operator_print_imag(sim->floquet_X[k_idx]);
        for (i = 0; i < qsn; i++)
        {
            for (j = 0; j < qsn; j++)
            {
                gsl_complex z;
                double delta, gamma;
                delta = ((qdpack_matrix_get(sim->floquet_quasienergies, i, 0)).dat[0])
                      - ((qdpack_matrix_get(sim->floquet_quasienergies, j, 0)).dat[0])
                      + 2 * 3.14159265358979323846 * k / T;
                do {(&z)->dat[0]=(delta); (&z)->dat[1]=(0);} while(0);
                qdpack_operator_set(sim->floquet_delta[k_idx], i, j, z);
                gamma = 2 * 3.14159265358979323846 * heaviside_function(delta) * s_func_cb(sim, delta)
                                 * gsl_complex_abs2(qdpack_operator_get(sim->floquet_X[k_idx], i, j));
                do {(&z)->dat[0]=(gamma); (&z)->dat[1]=(0);} while(0);
                qdpack_operator_set(sim->floquet_gamma[k_idx], i, j, z);
            }
        }
        printf("%s: floquet_delta[%d]\n", __PRETTY_FUNCTION__, k_idx);
        qdpack_operator_print_real(sim->floquet_delta[k_idx]);
        printf("%s: floquet_gamma[%d]\n", __PRETTY_FUNCTION__, k_idx);
        qdpack_operator_print_real(sim->floquet_gamma[k_idx]);
    }
    for (i = 0; i < qsn; i++)
    {
        for (j = 0; j < qsn; j++)
        {
            double a_ij_sum = 0.0;
            for (k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++)
            {
                double g1, g2, d;
                int k1_idx, k2_idx;
                k1_idx = k + sim->floquet_k_max;
                k2_idx = - k + sim->floquet_k_max;
                d = fabs(((qdpack_operator_get(sim->floquet_delta[k1_idx], i, j)).dat[0]));
                g1 = ((qdpack_operator_get(sim->floquet_gamma[k1_idx], i, j)).dat[0]);
                g2 = ((qdpack_operator_get(sim->floquet_gamma[k2_idx], j, i)).dat[0]);
                a_ij_sum += g1 + boltzmann_factor(d, sim->T_w) * (g1 + g2);
            }
            do {(&z)->dat[0]=(a_ij_sum); (&z)->dat[1]=(0);} while(0);
            qdpack_operator_set(sim->floquet_A, i, j, z);
        }
    }
    printf("%s: floquet_A\n", __PRETTY_FUNCTION__);
    qdpack_operator_print_real(sim->floquet_A);
    qdpack_matrix_set(sim->floquet_Gamma, 0, 1,
                      gsl_complex_rect(((qdpack_operator_get(sim->floquet_A, 0, 1)).dat[0]), 0.0));
    qdpack_matrix_set(sim->floquet_Gamma, 1, 0,
                      gsl_complex_rect(((qdpack_operator_get(sim->floquet_A, 1, 0)).dat[0]), 0.0));
    qdpack_matrix_set(sim->floquet_Gamma, 0, 0,
                      gsl_complex_rect(
                        (((qdpack_operator_get(sim->floquet_A, 0, 0)).dat[0]) +
                         ((qdpack_operator_get(sim->floquet_A, 0, 1)).dat[0]) +
                         ((qdpack_operator_get(sim->floquet_A, 1, 0)).dat[0]) +
                         ((qdpack_operator_get(sim->floquet_A, 1, 1)).dat[0])) / 2.0, 0.0));
    qdpack_matrix_set(sim->floquet_Gamma, 1, 1, qdpack_matrix_get(sim->floquet_Gamma, 0, 0));
    if (sim->floquet_ss_prob == ((void *)0))
    {
        sim->floquet_ss_prob = qdpack_matrix_alloc(qsn, 0);
    }
    qdpack_matrix_set(sim->floquet_ss_prob, 0, 0,
                      gsl_complex_rect(
                          ((qdpack_matrix_get(sim->floquet_Gamma, 1, 0)).dat[0]) /
                          (((qdpack_matrix_get(sim->floquet_Gamma, 0, 1)).dat[0]) +
                           ((qdpack_matrix_get(sim->floquet_Gamma, 1, 0)).dat[0])), 0.0));
    qdpack_matrix_set(sim->floquet_ss_prob, 1, 0,
                      gsl_complex_rect(
                          ((qdpack_matrix_get(sim->floquet_Gamma, 0, 1)).dat[0]) /
                          (((qdpack_matrix_get(sim->floquet_Gamma, 0, 1)).dat[0]) +
                           ((qdpack_matrix_get(sim->floquet_Gamma, 1, 0)).dat[0])), 0.0));
    for (k_idx = 0, k = - sim->floquet_k_max; k < sim->floquet_k_max; k++, k_idx++)
    {
        qdpack_operator_free(sim->floquet_X[k_idx]);
        qdpack_operator_free(sim->floquet_delta[k_idx]);
        qdpack_operator_free(sim->floquet_gamma[k_idx]);
    }
    qdpack_operator_free(w1);
    qdpack_operator_free(w2);
    qdpack_operator_free(sim->floquet_A);
    free(sim->floquet_X);
    free(sim->floquet_delta);
    free(sim->floquet_gamma);
    return 0;
}

canonical_original_code:
int qdpack_floquet_dissipation_rates(struct _qdpack_simulation *sim, struct qdpack_hilbert_space_t *qs, int (*floquet_modes_store_cb)(struct _qdpack_simulation *sp, struct qdpack_hilbert_space_t *qs, struct _qdpack_operator_t *rho_t, double t), double (*s_func_cb)(struct _qdpack_simulation *sp, double f))
{
    struct gsl_complex z;
    double t; double delta_t; double T;
    int i; int i_limit = 100; int j; int qsn = qdpack_hilbert_space_nstates(qs);
    struct _qdpack_operator_t *w1; struct _qdpack_operator_t *w2;
    int k; int k_idx;
    w1 = qdpack_operator_alloc(qs);
    w2 = qdpack_operator_alloc(qs);
    sim->floquet_k_max = 1;
    sim->floquet_gamma = malloc(sizeof(struct _qdpack_operator_t *) * (2 * sim->floquet_k_max + 1));
    sim->floquet_delta = malloc(sizeof(struct _qdpack_operator_t *) * (2 * sim->floquet_k_max + 1));
    sim->floquet_X = malloc(sizeof(struct _qdpack_operator_t *) * (2 * sim->floquet_k_max + 1));
    if (sim->floquet_X == ((void *)0))
    {
        fprintf(stderr, "%s: Failed to allocate floquet dissipation matrices.\n", __PRETTY_FUNCTION__);
        return -1;
    }
    for (k_idx = 0, k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++, k_idx++)
    {
        sim->floquet_X[k_idx] = qdpack_operator_alloc(qs);
        qdpack_operator_set_zero(sim->floquet_X[k_idx]);
        sim->floquet_delta[k_idx] = qdpack_operator_alloc(qs);
        qdpack_operator_set_zero(sim->floquet_delta[k_idx]);
        sim->floquet_gamma[k_idx] = qdpack_operator_alloc(qs);
        qdpack_operator_set_zero(sim->floquet_gamma[k_idx]);
    }
    sim->floquet_modes_t = qdpack_operator_alloc(qs);
    sim->floquet_states_t = qdpack_operator_alloc(qs);
    sim->floquet_A = qdpack_operator_alloc(qs);
    if (sim->floquet_Gamma == ((void *)0))
    {
        sim->floquet_Gamma = qdpack_matrix_alloc(qsn, qsn);
    }
    floquet_modes_store_cb(sim, qs, sim->floquet_modes, 0.0);
    T = 2*3.14159265358979323846/sim->h_td_w;
    delta_t = T / i_limit;
    for (t = delta_t; t <= T+delta_t/2.0; t += delta_t)
    {
        qdpack_floquet_modes_t(sim, qs, t);
        floquet_modes_store_cb(sim, qs, sim->floquet_modes_t, t);
        for (k_idx = 0, k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++, k_idx++)
        {
            qdpack_operator_blas_zgemm(CblasConjTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), sim->floquet_modes_t, sim->floquet_operator, (gsl_complex_rect(0.0,0.0)), w1);
            qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), w1 , sim->floquet_modes_t, (gsl_complex_rect(0.0,0.0)), w2);
            z = gsl_complex_mul_real(gsl_complex_polar(1.0, - 2 * 3.14159265358979323846 * k * t / T), delta_t / T);
            qdpack_operator_scale(w2, z);
            qdpack_operator_add(sim->floquet_X[k_idx], w2);
        }
    }
    for (k_idx = 0, k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++, k_idx++)
    {
        printf("%s: floquet_X[%d]\n", __PRETTY_FUNCTION__, k_idx);
        qdpack_operator_print_real(sim->floquet_X[k_idx]);
        qdpack_operator_print_imag(sim->floquet_X[k_idx]);
        for (i = 0; i < qsn; i++)
        {
            for (j = 0; j < qsn; j++)
            {
                struct gsl_complex z;
                double delta; double gamma;
                delta = ((qdpack_matrix_get(sim->floquet_quasienergies, i, 0)).dat[0])
                      - ((qdpack_matrix_get(sim->floquet_quasienergies, j, 0)).dat[0])
                      + 2 * 3.14159265358979323846 * k / T;
                do {(&z)->dat[0]=(delta); (&z)->dat[1]=(0);} while(0);
                qdpack_operator_set(sim->floquet_delta[k_idx], i, j, z);
                gamma = 2 * 3.14159265358979323846 * heaviside_function(delta) * s_func_cb(sim, delta)
                                 * gsl_complex_abs2(qdpack_operator_get(sim->floquet_X[k_idx], i, j));
                do {(&z)->dat[0]=(gamma); (&z)->dat[1]=(0);} while(0);
                qdpack_operator_set(sim->floquet_gamma[k_idx], i, j, z);
            }
        }
        printf("%s: floquet_delta[%d]\n", __PRETTY_FUNCTION__, k_idx);
        qdpack_operator_print_real(sim->floquet_delta[k_idx]);
        printf("%s: floquet_gamma[%d]\n", __PRETTY_FUNCTION__, k_idx);
        qdpack_operator_print_real(sim->floquet_gamma[k_idx]);
    }
    for (i = 0; i < qsn; i++)
    {
        for (j = 0; j < qsn; j++)
        {
            double a_ij_sum = 0.0;
            for (k = - sim->floquet_k_max; k <= sim->floquet_k_max; k++)
            {
                double g1; double g2; double d;
                int k1_idx; int k2_idx;
                k1_idx = k + sim->floquet_k_max;
                k2_idx = - k + sim->floquet_k_max;
                d = fabs(((qdpack_operator_get(sim->floquet_delta[k1_idx], i, j)).dat[0]));
                g1 = ((qdpack_operator_get(sim->floquet_gamma[k1_idx], i, j)).dat[0]);
                g2 = ((qdpack_operator_get(sim->floquet_gamma[k2_idx], j, i)).dat[0]);
                a_ij_sum += g1 + boltzmann_factor(d, sim->T_w) * (g1 + g2);
            }
            do {(&z)->dat[0]=(a_ij_sum); (&z)->dat[1]=(0);} while(0);
            qdpack_operator_set(sim->floquet_A, i, j, z);
        }
    }
    printf("%s: floquet_A\n", __PRETTY_FUNCTION__);
    qdpack_operator_print_real(sim->floquet_A);
    qdpack_matrix_set(sim->floquet_Gamma, 0, 1,
                      gsl_complex_rect(((qdpack_operator_get(sim->floquet_A, 0, 1)).dat[0]), 0.0));
    qdpack_matrix_set(sim->floquet_Gamma, 1, 0,
                      gsl_complex_rect(((qdpack_operator_get(sim->floquet_A, 1, 0)).dat[0]), 0.0));
    qdpack_matrix_set(sim->floquet_Gamma, 0, 0,
                      gsl_complex_rect(
                        (((qdpack_operator_get(sim->floquet_A, 0, 0)).dat[0]) +
                         ((qdpack_operator_get(sim->floquet_A, 0, 1)).dat[0]) +
                         ((qdpack_operator_get(sim->floquet_A, 1, 0)).dat[0]) +
                         ((qdpack_operator_get(sim->floquet_A, 1, 1)).dat[0])) / 2.0, 0.0));
    qdpack_matrix_set(sim->floquet_Gamma, 1, 1, qdpack_matrix_get(sim->floquet_Gamma, 0, 0));
    if (sim->floquet_ss_prob == ((void *)0))
    {
        sim->floquet_ss_prob = qdpack_matrix_alloc(qsn, 0);
    }
    qdpack_matrix_set(sim->floquet_ss_prob, 0, 0,
                      gsl_complex_rect(
                          ((qdpack_matrix_get(sim->floquet_Gamma, 1, 0)).dat[0]) /
                          (((qdpack_matrix_get(sim->floquet_Gamma, 0, 1)).dat[0]) +
                           ((qdpack_matrix_get(sim->floquet_Gamma, 1, 0)).dat[0])), 0.0));
    qdpack_matrix_set(sim->floquet_ss_prob, 1, 0,
                      gsl_complex_rect(
                          ((qdpack_matrix_get(sim->floquet_Gamma, 0, 1)).dat[0]) /
                          (((qdpack_matrix_get(sim->floquet_Gamma, 0, 1)).dat[0]) +
                           ((qdpack_matrix_get(sim->floquet_Gamma, 1, 0)).dat[0])), 0.0));
    for (k_idx = 0, k = - sim->floquet_k_max; k < sim->floquet_k_max; k++, k_idx++)
    {
        qdpack_operator_free(sim->floquet_X[k_idx]);
        qdpack_operator_free(sim->floquet_delta[k_idx]);
        qdpack_operator_free(sim->floquet_gamma[k_idx]);
    }
    qdpack_operator_free(w1);
    qdpack_operator_free(w2);
    qdpack_operator_free(sim->floquet_A);
    free(sim->floquet_X);
    free(sim->floquet_delta);
    free(sim->floquet_gamma);
    return 0;
}

variable_types: int (*)(struct _qdpack_simulation *sp, struct qdpack_hilbert_space_t *qs, struct _qdpack_operator_t *rho_t, double t)

variable_types: double (*)(struct _qdpack_simulation *sp, double f)

variable_types:
z
struct gsl_complex { double dat[2]; }

variable_types: T double

variable_types: k int

variable_types: delta double

variable_types: gamma double

variable_types: d double


original_code:
qdpack_operator_t *
basis_transform(qdpack_simulation_t *sim, qdpack_hilbert_space_t *qs, int n, ho_w_func_t ho_w_cb, double t)
{
    gsl_complex z;
    qdpack_operator_t *S;
    int i, j, ns, mm, ll;
    double r, w0, wt, sij;
    qdpack_operator_list_t op_list;
    ns = qs->nstates[n];
    if ((S = qdpack_operator_alloc(qs)) == ((void *)0))
    {
        fprintf(stderr, "basis_transform_ho: Failed to allocate basis tranformation matrix (Out of memory?).\n");
    }
    qdpack_operator_set_zero(S);
    w0 = ho_w_cb(0.0, sim);
    wt = ho_w_cb(t, sim);
    r = 2 * sqrt(w0 * wt) / (wt + w0);
    if (r > 1.0)
        r = 1.0;
    else if (r < -1.0)
        r = -1.0;
    for (j = 0; j < ns; j++)
    {
        for (i = j; i < ns; i++)
        {
            if ((i == j) || ((i+j)%2 == 0))
            {
                ll = (i+j)/2;
                mm = abs(i-j)/2;
                sij = fact_sqrt(((i<j)?i:j)) / fact_sqrt(((i<j)?j:i)) * sqrt(r) * gsl_sf_legendre_Plm(ll, mm, r);
                do {(&z)->dat[0]=(sij); (&z)->dat[1]=(0.0);} while(0);
                qdpack_operator_set(S, i, j, z);
                if (i != j)
                {
                    sij = pow(-1.0, (j-i)/2.0) * sij;
                    do {(&z)->dat[0]=(sij); (&z)->dat[1]=(0.0);} while(0);
                    qdpack_operator_set(S, j, i, z);
                }
            }
        }
    }
    if (qs->nsubsys > 1)
    {
        qdpack_operator_list_init(&op_list);
        for (i = 0; i < qs->nsubsys; i++)
        {
            qdpack_operator_t *M;
            qdpack_hilbert_space_t *qs_sub = qdpack_hilbert_space_new();
            qdpack_hilbert_space_add(qs_sub, qs->nstates[i]);
            M = qdpack_operator_alloc(qs_sub);
            if (i != n)
            {
                operator_unit(M, qs, i);
            }
            else
            {
                M = S;
            }
            qdpack_operator_list_append(&op_list, M);
        }
        qdpack_operator_tensor(S, qs, &op_list);
        qdpack_operator_list_free_entries(&op_list);
    }
    return S;
}

canonical_original_code:
struct _qdpack_operator_t *basis_transform(struct _qdpack_simulation *sim, struct qdpack_hilbert_space_t *qs, int n, double (*ho_w_cb)(double, void *), double t)
{
    struct gsl_complex z;
    struct _qdpack_operator_t *S;
    int i; int j; int ns; int mm; int ll;
    double r; double w0; double wt; double sij;
    struct qdpack_operator_list_t op_list;
    ns = qs->nstates[n];
    if ((S = qdpack_operator_alloc(qs)) == ((void *)0))
    {
        fprintf(stderr, "basis_transform_ho: Failed to allocate basis tranformation matrix (Out of memory?).\n");
    }
    qdpack_operator_set_zero(S);
    w0 = ho_w_cb(0.0, sim);
    wt = ho_w_cb(t, sim);
    r = 2 * sqrt(w0 * wt) / (wt + w0);
    if (r > 1.0)
        r = 1.0;
    else if (r < -1.0)
        r = -1.0;
    for (j = 0; j < ns; j++)
    {
        for (i = j; i < ns; i++)
        {
            if ((i == j) || ((i+j)%2 == 0))
            {
                ll = (i+j)/2;
                mm = abs(i-j)/2;
                sij = fact_sqrt(((i<j)?i:j)) / fact_sqrt(((i<j)?j:i)) * sqrt(r) * gsl_sf_legendre_Plm(ll, mm, r);
                do {(&z)->dat[0]=(sij); (&z)->dat[1]=(0.0);} while(0);
                qdpack_operator_set(S, i, j, z);
                if (i != j)
                {
                    sij = pow(-1.0, (j-i)/2.0) * sij;
                    do {(&z)->dat[0]=(sij); (&z)->dat[1]=(0.0);} while(0);
                    qdpack_operator_set(S, j, i, z);
                }
            }
        }
    }
    if (qs->nsubsys > 1)
    {
        qdpack_operator_list_init(&op_list);
        for (i = 0; i < qs->nsubsys; i++)
        {
            struct _qdpack_operator_t *M;
            struct qdpack_hilbert_space_t *qs_sub = qdpack_hilbert_space_new();
            qdpack_hilbert_space_add(qs_sub, qs->nstates[i]);
            M = qdpack_operator_alloc(qs_sub);
            if (i != n)
            {
                operator_unit(M, qs, i);
            }
            else
            {
                M = S;
            }
            qdpack_operator_list_append(&op_list, M);
        }
        qdpack_operator_tensor(S, qs, &op_list);
        qdpack_operator_list_free_entries(&op_list);
    }
    return S;
}

variable_types: struct _qdpack_simulation { struct qdpack_hilbert_space_t *qs; double Ti; double Tf; double dT; char simsig[1024]; double T_w; struct _qdpack_operator_t *H0; struct _qdpack_operator_t *H1; struct _qdpack_operator_t *H_t; struct _qdpack_operator_t *H_a; struct _qdpack_operator_t *H_b; struct _qdpack_operator_t *do_a[20]; struct _qdpack_operator_t *do_ad[20]; struct _qdpack_operator_t *do_aad[20]; struct _qdpack_operator_t *do_ada[20]; struct _qdpack_operator_t *do_a_eb[20]; struct _qdpack_operator_t *do_ad_eb[20]; struct _qdpack_operator_t *do_aad_eb[20]; struct _qdpack_operator_t *do_ada_eb[20]; struct _qdpack_operator_t *ws1; struct _qdpack_operator_t *ws2; struct _qdpack_operator_t *op_tmp; double do_g1[20]; double do_g2[20]; int do_n; void *parameters; struct _qdpack_operator_t *dSSt; double (*ho_w_func)(double, void *); int cont_basis_transform; int option_me_eb; int option_adaptive; double g0; struct _qdpack_operator_t *N_op[10]; struct _qdpack_operator_t *a; struct _qdpack_operator_t *n_op; struct _qdpack_state { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix *data; } *psi_f; double h_td_A; double h_td_w; int (*H0_func)(struct _qdpack_simulation *, struct qdpack_hilbert_space_t *, struct _qdpack_operator_t *); int (*H1_func)(struct _qdpack_simulation *, struct qdpack_hilbert_space_t *, struct _qdpack_operator_t *); int (*Ht_func)(struct _qdpack_simulation *, struct _qdpack_operator_t *, struct _qdpack_operator_t *, struct _qdpack_operator_t *, double); struct _qdpack_operator_t *evec; struct _qdpack_matrix *eval; struct _qdpack_operator_t *floquet_states_t; struct _qdpack_operator_t *floquet_modes; struct _qdpack_matrix *floquet_quasienergies; struct _qdpack_operator_t **floquet_X; struct _qdpack_matrix *floquet_Gamma; struct _qdpack_matrix *floquet_ss_prob; int floquet_k_max; } *

variable_types: n int

variable_types: ns int

variable_types: mm int

variable_types: ll int

variable_types: r double

variable_types: wt double

variable_types: sij double

variable_types: struct qdpack_operator_list_t { struct _qdpack_operator_t { int n; struct qdpack_hilbert_space_t { int nsubsys; int *nstates; } *qs; struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *data; } *dm[20]; int len; }

user_defined_types: struct qdpack_operator_list_t { struct _qdpack_operator_t *dm[20]; int len; }


decompiled_code:
__int64 __fastcall basis_transform_derivative(__int64 a1, __int64 a2, int a3, double (__fastcall *a4)(__int64, double), double a5, double a6, double a7, double a8)
{
__int64 v10; // r13
__int64 v11; // r12
double v12; // r14
__int64 v13; // xmm0_8

v10 = basis_transform(a1, a2, a3, a4, a5 - 0.000001, a6, a7, a8);
v11 = basis_transform(a1, a2, a3, a4, a5 + 0.000001, a6, a7, a5 + 0.000001);
qdpack_operator_sub(v11, v10);
v12 = a4(a1, a5);
if ( a4(a1, 0.0) <= v12 )
v13 = LC8;
else
v13 = LC7;
qdpack_operator_scale(v11, *(double *)&v13, 0.0);
qdpack_operator_free(v10);
return v11;
}

original_code:
qdpack_operator_t *
basis_transform_derivative(qdpack_simulation_t *sim, qdpack_hilbert_space_t *qs, int n, ho_w_func_t ho_w_cb, double t)
{
    double h;
    gsl_complex z;
    qdpack_operator_t *S1, *S2;
    h = 1e-6;
    S1 = basis_transform(sim, qs, n, ho_w_cb, t-h);
    S2 = basis_transform(sim, qs, n, ho_w_cb, t+h);
    qdpack_operator_sub(S2, S1);
    if (ho_w_cb(t, sim) < ho_w_cb(0.0, sim))
    {
        do {(&z)->dat[0]=(-1.0 / (2.0 * h)); (&z)->dat[1]=(0.0);} while(0);
    }
    else
    {
        do {(&z)->dat[0]=(1.0 / (2.0 * h)); (&z)->dat[1]=(0.0);} while(0);
    }
    qdpack_operator_scale(S2, z);
    qdpack_operator_free(S1);
    return S2;
}

canonical_original_code:
struct _qdpack_operator_t *basis_transform_derivative(struct _qdpack_simulation *sim, struct qdpack_hilbert_space_t *qs, int n, double (*ho_w_cb)(double, void *), double t)
{
    double h;
    struct gsl_complex z;
    struct _qdpack_operator_t *S1; struct _qdpack_operator_t *S2;
    h = 1e-6;
    S1 = basis_transform(sim, qs, n, ho_w_cb, t-h);
    S2 = basis_transform(sim, qs, n, ho_w_cb, t+h);
    qdpack_operator_sub(S2, S1);
    if (ho_w_cb(t, sim) < ho_w_cb(0.0, sim))
    {
        do {(&z)->dat[0]=(-1.0 / (2.0 * h)); (&z)->dat[1]=(0.0);} while(0);
    }
    else
    {
        do {(&z)->dat[0]=(1.0 / (2.0 * h)); (&z)->dat[1]=(0.0);} while(0);
    }
    qdpack_operator_scale(S2, z);
    qdpack_operator_free(S1);
    return S2;
}

variable_types:
qs
struct qdpack_hilbert_space_t *

variable_types: t double

variable_types: h double

variable_types: struct _qdpack_operator_t { int n; struct qdpack_hilbert_space_t { int nsubsys; int *nstates; } *qs; struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *data; } *


decompiled_code:
double __fastcall qdpack_entanglement_neumann_entropy(int *a1)
{
__int64 v1; // r13
__int64 v2; // r12
__int64 v3; // r14
__int64 v4; // rbx
double v5; // r15
double v7; // [rsp+8h] [rbp-40h]

v1 = qdpack_matrix_alloc(a1[1], *a1);
qdpack_matrix_memcpy(v1, *((_QWORD *)a1 + 2));
v2 = qdpack_matrix_alloc(a1[1], 1LL);
v3 = qdpack_matrix_alloc(a1[1], *a1);
qdpack_matrix_eigen_hermv(v1, v2, v3, 0LL);
if ( a1[1] <= 0 )
{
v5 = *(double *)&LC0;
}
else
{
v4 = 0LL;
v5 = *(double *)&LC0;
do
{
v7 = qdpack_matrix_get(v2, v4, 0LL);
v5 = v5 - log(v7) * v7;
++v4;
}
while ( a1[1] > (int)v4 );
}
qdpack_matrix_free(v1);
qdpack_matrix_free(v2);
qdpack_matrix_free(v3);
return v5;
}

canonical_decompiled_code:
double __fastcall func0(int *a1)
{
__int64 v1; // r13
__int64 v2; // r12
__int64 v3; // r14
__int64 v4; // rbx
double v5; // r15
double v7; // [rsp+8h] [rbp-40h]

v1 = qdpack_matrix_alloc(a1[1], *a1);
qdpack_matrix_memcpy(v1, *((_QWORD *)a1 + 2));
v2 = qdpack_matrix_alloc(a1[1], 1LL);
v3 = qdpack_matrix_alloc(a1[1], *a1);
qdpack_matrix_eigen_hermv(v1, v2, v3, 0LL);
if ( a1[1] <= 0 )
{
v5 = *(double *)&LC0;
}
else
{
v4 = 0LL;
v5 = *(double *)&LC0;
do
{
v7 = qdpack_matrix_get(v2, v4, 0LL);
v5 = v5 - log(v7) * v7;
++v4;
}
while ( a1[1] > (int)v4 );
}
qdpack_matrix_free(v1);
qdpack_matrix_free(v2);
qdpack_matrix_free(v3);
return v5;
}

original_code:
double
qdpack_entanglement_neumann_entropy(qdpack_operator_t *rho)
{
    int i;
    double e, ne = 0.0;
    qdpack_matrix_t *rho_w;
    qdpack_matrix_t *eval, *evec;
    rho_w = qdpack_matrix_alloc(rho->m, rho->n);
    qdpack_matrix_memcpy(rho_w, rho->data);
    eval = qdpack_matrix_alloc(rho->m, 1);
    evec = qdpack_matrix_alloc(rho->m, rho->n);
    qdpack_matrix_eigen_hermv(rho_w, eval, evec, 0);
    for (i = 0; i < rho->m; i++)
    {
        e = ((qdpack_matrix_get(eval, i, 0)).dat[0]);
        ne += - e * log(e);
    }
    qdpack_matrix_free(rho_w);
    qdpack_matrix_free(eval);
    qdpack_matrix_free(evec);
    return ne;
}

canonical_original_code:
double qdpack_entanglement_neumann_entropy(struct _qdpack_operator_t *rho)
{
    int i;
    double e; double ne = 0.0;
    struct _qdpack_matrix *rho_w;
    struct _qdpack_matrix *eval; struct _qdpack_matrix *evec;
    rho_w = qdpack_matrix_alloc(rho->m, rho->n);
    qdpack_matrix_memcpy(rho_w, rho->data);
    eval = qdpack_matrix_alloc(rho->m, 1);
    evec = qdpack_matrix_alloc(rho->m, rho->n);
    qdpack_matrix_eigen_hermv(rho_w, eval, evec, 0);
    for (i = 0; i < rho->m; i++)
    {
        e = ((qdpack_matrix_get(eval, i, 0)).dat[0]);
        ne += - e * log(e);
    }
    qdpack_matrix_free(rho_w);
    qdpack_matrix_free(eval);
    qdpack_matrix_free(evec);
    return ne;
}

variable_types: i int

variable_types: e double

variable_types: ne double

variable_types: struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *

return_type: double


decompiled_code:
double __fastcall qdpack_entanglement_log_neg(__int64 a1)
{
__int64 v1; // rbp
__int64 v2; // r12
__int64 v3; // r13
__int64 v4; // rbx
double v5; // r14
double v6; // xmm0_8
double v7; // xmm0_8
double v8; // rbx

v1 = qdpack_matrix_alloc(4LL, 4LL);
qdpack_matrix_memcpy(v1, *(_QWORD *)(a1 + 16));
qdpack_operator_get(a1, 2LL, 0LL);
qdpack_matrix_set(v1, 0LL, 2LL);
qdpack_operator_get(a1, 2LL, 1LL);
qdpack_matrix_set(v1, 0LL, 3LL);
qdpack_operator_get(a1, 3LL, 0LL);
qdpack_matrix_set(v1, 1LL, 2LL);
qdpack_operator_get(a1, 3LL, 1LL);
qdpack_matrix_set(v1, 1LL, 3LL);
qdpack_operator_get(a1, 0LL, 2LL);
qdpack_matrix_set(v1, 2LL, 0LL);
qdpack_operator_get(a1, 0LL, 3LL);
qdpack_matrix_set(v1, 2LL, 1LL);
qdpack_operator_get(a1, 1LL, 2LL);
qdpack_matrix_set(v1, 3LL, 0LL);
qdpack_operator_get(a1, 1LL, 3LL);
qdpack_matrix_set(v1, 3LL, 1LL);
v2 = qdpack_matrix_alloc(4LL, 1LL);
v3 = qdpack_matrix_alloc(4LL, 4LL);
qdpack_matrix_eigen_hermv(v1, v2, v3, 0LL);
v4 = 0LL;
v5 = 0.0;
do
{
v6 = qdpack_matrix_get(v2, v4, 0LL);
v7 = gsl_complex_abs(v6) + v5;
v5 = v7;
++v4;
}
while ( v4 != 4 );
v8 = log2(v7);
qdpack_matrix_free(v1);
qdpack_matrix_free(v2);
qdpack_matrix_free(v3);
return v8;
}

canonical_decompiled_code:
double __fastcall func1(__int64 a1)
{
__int64 v1; // rbp
__int64 v2; // r12
__int64 v3; // r13
__int64 v4; // rbx
double v5; // r14
double v6; // xmm0_8
double v7; // xmm0_8
double v8; // rbx

v1 = qdpack_matrix_alloc(4LL, 4LL);
qdpack_matrix_memcpy(v1, *(_QWORD *)(a1 + 16));
qdpack_operator_get(a1, 2LL, 0LL);
qdpack_matrix_set(v1, 0LL, 2LL);
qdpack_operator_get(a1, 2LL, 1LL);
qdpack_matrix_set(v1, 0LL, 3LL);
qdpack_operator_get(a1, 3LL, 0LL);
qdpack_matrix_set(v1, 1LL, 2LL);
qdpack_operator_get(a1, 3LL, 1LL);
qdpack_matrix_set(v1, 1LL, 3LL);
qdpack_operator_get(a1, 0LL, 2LL);
qdpack_matrix_set(v1, 2LL, 0LL);
qdpack_operator_get(a1, 0LL, 3LL);
qdpack_matrix_set(v1, 2LL, 1LL);
qdpack_operator_get(a1, 1LL, 2LL);
qdpack_matrix_set(v1, 3LL, 0LL);
qdpack_operator_get(a1, 1LL, 3LL);
qdpack_matrix_set(v1, 3LL, 1LL);
v2 = qdpack_matrix_alloc(4LL, 1LL);
v3 = qdpack_matrix_alloc(4LL, 4LL);
qdpack_matrix_eigen_hermv(v1, v2, v3, 0LL);
v4 = 0LL;
v5 = 0.0;
do
{
v6 = qdpack_matrix_get(v2, v4, 0LL);
v7 = gsl_complex_abs(v6) + v5;
v5 = v7;
++v4;
}
while ( v4 != 4 );
v8 = log2(v7);
qdpack_matrix_free(v1);
qdpack_matrix_free(v2);
qdpack_matrix_free(v3);
return v8;
}

original_code:
double
qdpack_entanglement_log_neg(qdpack_operator_t *rho)
{
    int i;
    double ln, sum_eval;
    qdpack_matrix_t *rho_pt;
    qdpack_matrix_t *eval, *evec;
    rho_pt = qdpack_matrix_alloc(4,4);
    qdpack_matrix_memcpy(rho_pt, rho->data);
    qdpack_matrix_set(rho_pt, 0, 2, qdpack_operator_get(rho, 2, 0));
    qdpack_matrix_set(rho_pt, 0, 3, qdpack_operator_get(rho, 2, 1));
    qdpack_matrix_set(rho_pt, 1, 2, qdpack_operator_get(rho, 3, 0));
    qdpack_matrix_set(rho_pt, 1, 3, qdpack_operator_get(rho, 3, 1));
    qdpack_matrix_set(rho_pt, 2, 0, qdpack_operator_get(rho, 0, 2));
    qdpack_matrix_set(rho_pt, 2, 1, qdpack_operator_get(rho, 0, 3));
    qdpack_matrix_set(rho_pt, 3, 0, qdpack_operator_get(rho, 1, 2));
    qdpack_matrix_set(rho_pt, 3, 1, qdpack_operator_get(rho, 1, 3));
    eval = qdpack_matrix_alloc(4, 1);
    evec = qdpack_matrix_alloc(4, 4);
    qdpack_matrix_eigen_hermv(rho_pt, eval, evec, 0);
    sum_eval = 0;
    for (i = 0; i < 4; i++)
    {
        sum_eval += gsl_complex_abs(qdpack_matrix_get(eval, i, 0));
    }
    ln = log2(sum_eval);
    qdpack_matrix_free(rho_pt);
    qdpack_matrix_free(eval);
    qdpack_matrix_free(evec);
    return ln;
}

canonical_original_code:
double qdpack_entanglement_log_neg(struct _qdpack_operator_t *rho)
{
    int i;
    double ln; double sum_eval;
    struct _qdpack_matrix *rho_pt;
    struct _qdpack_matrix *eval; struct _qdpack_matrix *evec;
    rho_pt = qdpack_matrix_alloc(4,4);
    qdpack_matrix_memcpy(rho_pt, rho->data);
    qdpack_matrix_set(rho_pt, 0, 2, qdpack_operator_get(rho, 2, 0));
    qdpack_matrix_set(rho_pt, 0, 3, qdpack_operator_get(rho, 2, 1));
    qdpack_matrix_set(rho_pt, 1, 2, qdpack_operator_get(rho, 3, 0));
    qdpack_matrix_set(rho_pt, 1, 3, qdpack_operator_get(rho, 3, 1));
    qdpack_matrix_set(rho_pt, 2, 0, qdpack_operator_get(rho, 0, 2));
    qdpack_matrix_set(rho_pt, 2, 1, qdpack_operator_get(rho, 0, 3));
    qdpack_matrix_set(rho_pt, 3, 0, qdpack_operator_get(rho, 1, 2));
    qdpack_matrix_set(rho_pt, 3, 1, qdpack_operator_get(rho, 1, 3));
    eval = qdpack_matrix_alloc(4, 1);
    evec = qdpack_matrix_alloc(4, 4);
    qdpack_matrix_eigen_hermv(rho_pt, eval, evec, 0);
    sum_eval = 0;
    for (i = 0; i < 4; i++)
    {
        sum_eval += gsl_complex_abs(qdpack_matrix_get(eval, i, 0));
    }
    ln = log2(sum_eval);
    qdpack_matrix_free(rho_pt);
    qdpack_matrix_free(eval);
    qdpack_matrix_free(evec);
    return ln;
}

variable_types: ln double


decompiled_code:
double __fastcall qdpack_entanglement_concurrence(__int64 a1, __int64 a2)
{
__int64 v2; // rbx
__int64 v3; // rbp
__int64 v4; // r12
__int64 v5; // r15
__int64 v6; // r12
double v7; // r15
double v8; // xmm0_8
double v9; // r15
double v10; // xmm0_8
double v11; // r15
double v12; // xmm0_8
double v13; // r14
double v14; // xmm0_8
__int64 i; // r13
double v16; // xmm0_8
__int64 j; // r14
__int64 k; // r13
double v19; // xmm0_8
double v20; // r14
double v21; // xmm0_8
double v22; // r13
double v23; // r13
double v24; // r13
double v25; // r13
__int64 v27; // [rsp+0h] [rbp-78h]
__int64 v28; // [rsp+8h] [rbp-70h]
__int64 v29; // [rsp+10h] [rbp-68h]
__int64 v30; // [rsp+18h] [rbp-60h]
double v31; // [rsp+28h] [rbp-50h]
double v32; // [rsp+38h] [rbp-40h]

v2 = qdpack_matrix_alloc(4LL, 1LL);
v3 = qdpack_matrix_alloc(4LL, 4LL);
v4 = qdpack_operator_alloc(a1);
v5 = qdpack_operator_alloc(a1);
v29 = v4;
operator_sigma_y(v4, a1, 0LL);
v30 = v5;
operator_sigma_y(v5, a1, 1LL);
v6 = qdpack_operator_alloc(a1);
v27 = qdpack_operator_alloc(a1);
v28 = qdpack_operator_alloc(a1);
v7 = gsl_complex_rect(0.0);
v8 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v29, v30, v6, v8, 0.0, v7, 0.0);
v9 = gsl_complex_rect(0.0);
v10 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(111LL, 111LL, a2, v6, v27, v10, 0.0, v9, 0.0);
v11 = gsl_complex_rect(0.0);
v12 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(113LL, 111LL, a2, v6, v28, v12, 0.0, v11, 0.0);
v13 = gsl_complex_rect(0.0);
v14 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v27, v28, v6, v14, 0.0, v13, 0.0);
puts("rr = ");
qdpack_operator_print(v6);
qdpack_matrix_eigen_zgeev(*(_QWORD *)(v6 + 16), v2, v3, 0LL);
puts("Eigenvalues = ");
for ( i = 0LL; i != 4; ++i )
{
v16 = qdpack_matrix_get(v2, i, 0LL);
printf("\t(%f, %f)", v16, 0.0);
}
putchar(10);
puts("evec = ");
for ( j = 0LL; j != 4; ++j )
{
for ( k = 0LL; k != 4; ++k )
{
v19 = qdpack_matrix_get(v3, j, k);
printf("\t(%f, %f)", v19, 0.0);
}
putchar(10);
}
qdpack_matrix_get(v2, 3LL, 0LL);
v20 = qdpack_matrix_get(v2, 3LL, 0LL);
qdpack_matrix_get(v2, 2LL, 0LL);
v31 = qdpack_matrix_get(v2, 2LL, 0LL);
qdpack_matrix_get(v2, 1LL, 0LL);
v32 = qdpack_matrix_get(v2, 1LL, 0LL);
qdpack_matrix_get(v2, 0LL, 0LL);
v21 = qdpack_matrix_get(v2, 0LL, 0LL);
printf("c: eigvals: (%f, %f), (%f, %f), (%f, %f), (%f, %f)\n", v21, 0.0, v32, 0.0, v31, 0.0, v20, 0.0);
v22 = qdpack_matrix_get(v2, 0LL, 0LL);
v23 = v22 - qdpack_matrix_get(v2, 1LL, 0LL);
v24 = v23 - qdpack_matrix_get(v2, 2LL, 0LL);
v25 = v24 - qdpack_matrix_get(v2, 3LL, 0LL);
if ( v25 < 0.0 )
v25 = 0.0;
qdpack_operator_free(v29);
qdpack_operator_free(v30);
qdpack_operator_free(v6);
qdpack_operator_free(v27);
qdpack_operator_free(v28);
qdpack_matrix_free(v3);
qdpack_matrix_free(v2);
return v25;
}

canonical_decompiled_code:
double __fastcall func2(__int64 a1, __int64 a2)
{
__int64 v2; // rbx
__int64 v3; // rbp
__int64 v4; // r12
__int64 v5; // r15
__int64 v6; // r12
double v7; // r15
double v8; // xmm0_8
double v9; // r15
double v10; // xmm0_8
double v11; // r15
double v12; // xmm0_8
double v13; // r14
double v14; // xmm0_8
__int64 i; // r13
double v16; // xmm0_8
__int64 j; // r14
__int64 k; // r13
double v19; // xmm0_8
double v20; // r14
double v21; // xmm0_8
double v22; // r13
double v23; // r13
double v24; // r13
double v25; // r13
__int64 v27; // [rsp+0h] [rbp-78h]
__int64 v28; // [rsp+8h] [rbp-70h]
__int64 v29; // [rsp+10h] [rbp-68h]
__int64 v30; // [rsp+18h] [rbp-60h]
double v31; // [rsp+28h] [rbp-50h]
double v32; // [rsp+38h] [rbp-40h]

v2 = qdpack_matrix_alloc(4LL, 1LL);
v3 = qdpack_matrix_alloc(4LL, 4LL);
v4 = qdpack_operator_alloc(a1);
v5 = qdpack_operator_alloc(a1);
v29 = v4;
operator_sigma_y(v4, a1, 0LL);
v30 = v5;
operator_sigma_y(v5, a1, 1LL);
v6 = qdpack_operator_alloc(a1);
v27 = qdpack_operator_alloc(a1);
v28 = qdpack_operator_alloc(a1);
v7 = gsl_complex_rect(0.0);
v8 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v29, v30, v6, v8, 0.0, v7, 0.0);
v9 = gsl_complex_rect(0.0);
v10 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(111LL, 111LL, a2, v6, v27, v10, 0.0, v9, 0.0);
v11 = gsl_complex_rect(0.0);
v12 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(113LL, 111LL, a2, v6, v28, v12, 0.0, v11, 0.0);
v13 = gsl_complex_rect(0.0);
v14 = gsl_complex_rect(1.0);
qdpack_operator_blas_zgemm(111LL, 111LL, v27, v28, v6, v14, 0.0, v13, 0.0);
puts("rr = ");
qdpack_operator_print(v6);
qdpack_matrix_eigen_zgeev(*(_QWORD *)(v6 + 16), v2, v3, 0LL);
puts("Eigenvalues = ");
for ( i = 0LL; i != 4; ++i )
{
v16 = qdpack_matrix_get(v2, i, 0LL);
printf("\t(%f, %f)", v16, 0.0);
}
putchar(10);
puts("evec = ");
for ( j = 0LL; j != 4; ++j )
{
for ( k = 0LL; k != 4; ++k )
{
v19 = qdpack_matrix_get(v3, j, k);
printf("\t(%f, %f)", v19, 0.0);
}
putchar(10);
}
qdpack_matrix_get(v2, 3LL, 0LL);
v20 = qdpack_matrix_get(v2, 3LL, 0LL);
qdpack_matrix_get(v2, 2LL, 0LL);
v31 = qdpack_matrix_get(v2, 2LL, 0LL);
qdpack_matrix_get(v2, 1LL, 0LL);
v32 = qdpack_matrix_get(v2, 1LL, 0LL);
qdpack_matrix_get(v2, 0LL, 0LL);
v21 = qdpack_matrix_get(v2, 0LL, 0LL);
printf("c: eigvals: (%f, %f), (%f, %f), (%f, %f), (%f, %f)\n", v21, 0.0, v32, 0.0, v31, 0.0, v20, 0.0);
v22 = qdpack_matrix_get(v2, 0LL, 0LL);
v23 = v22 - qdpack_matrix_get(v2, 1LL, 0LL);
v24 = v23 - qdpack_matrix_get(v2, 2LL, 0LL);
v25 = v24 - qdpack_matrix_get(v2, 3LL, 0LL);
if ( v25 < 0.0 )
v25 = 0.0;
qdpack_operator_free(v29);
qdpack_operator_free(v30);
qdpack_operator_free(v6);
qdpack_operator_free(v27);
qdpack_operator_free(v28);
qdpack_matrix_free(v3);
qdpack_matrix_free(v2);
return v25;
}

original_code:
double
qdpack_entanglement_concurrence(qdpack_hilbert_space_t *qs, qdpack_operator_t *rho)
{
    double c = 0.0;
    qdpack_operator_t *y1, *y2, *y, *ws1, *ws2;
    qdpack_matrix_t *eval = qdpack_matrix_alloc(4, 1);
    qdpack_matrix_t *evec = qdpack_matrix_alloc(4, 4);
    y1 = qdpack_operator_alloc(qs);
    y2 = qdpack_operator_alloc(qs);
    operator_sigma_y(y1, qs, 0);
    operator_sigma_y(y2, qs, 1);
    y = qdpack_operator_alloc(qs);
    ws1 = qdpack_operator_alloc(qs);
    ws2 = qdpack_operator_alloc(qs);
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), y1, y2, (gsl_complex_rect(0.0,0.0)), y);
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), rho, y, (gsl_complex_rect(0.0,0.0)), ws1);
    qdpack_operator_blas_zgemm(CblasConjTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), rho, y, (gsl_complex_rect(0.0,0.0)), ws2);
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), ws1, ws2, (gsl_complex_rect(0.0,0.0)), y);
    printf("rr = \n"); qdpack_operator_print(y);
    qdpack_matrix_eigen_zgeev(y->data, eval, evec, 0);
    {
        int i, j, n = 4;
        printf("Eigenvalues = \n");
        for (i = 0; i < n; i++)
        {
                gsl_complex z;
                z = qdpack_matrix_get(eval, i, 0);
                printf("\t(%f, %f)", ((z).dat[0]), ((z).dat[1]));
        }
        printf("\n");
        printf("evec = \n");
        for (i = 0; i < n; i++)
        {
                for (j = 0; j < n; j++)
                {
                        gsl_complex z;
                        z = qdpack_matrix_get(evec, i, j);
                        printf("\t(%f, %f)", ((z).dat[0]), ((z).dat[1]));
                }
                printf("\n");
        }
    }
    printf("c: eigvals: (%f, %f), (%f, %f), (%f, %f), (%f, %f)\n",
           ((qdpack_matrix_get(eval, 0, 0)).dat[0]), ((qdpack_matrix_get(eval, 0, 0)).dat[1]),
           ((qdpack_matrix_get(eval, 1, 0)).dat[0]), ((qdpack_matrix_get(eval, 1, 0)).dat[1]),
           ((qdpack_matrix_get(eval, 2, 0)).dat[0]), ((qdpack_matrix_get(eval, 2, 0)).dat[1]),
           ((qdpack_matrix_get(eval, 3, 0)).dat[0]), ((qdpack_matrix_get(eval, 3, 0)).dat[1]));
    c = + ((qdpack_matrix_get(eval, 0, 0)).dat[0]) - ((qdpack_matrix_get(eval, 1, 0)).dat[0])
        - ((qdpack_matrix_get(eval, 2, 0)).dat[0]) - ((qdpack_matrix_get(eval, 3, 0)).dat[0]);
    c = (c < 0.0) ? 0.0 : c;
    qdpack_operator_free(y1);
    qdpack_operator_free(y2);
    qdpack_operator_free(y);
    qdpack_operator_free(ws1);
    qdpack_operator_free(ws2);
    qdpack_matrix_free(evec);
    qdpack_matrix_free(eval);
    return c;
}

canonical_original_code:
double qdpack_entanglement_concurrence(struct qdpack_hilbert_space_t *qs, struct _qdpack_operator_t *rho)
{
    double c = 0.0;
    struct _qdpack_operator_t *y1; struct _qdpack_operator_t *y2; struct _qdpack_operator_t *y; struct _qdpack_operator_t *ws1; struct _qdpack_operator_t *ws2;
    struct _qdpack_matrix *eval = qdpack_matrix_alloc(4, 1);
    struct _qdpack_matrix *evec = qdpack_matrix_alloc(4, 4);
    y1 = qdpack_operator_alloc(qs);
    y2 = qdpack_operator_alloc(qs);
    operator_sigma_y(y1, qs, 0);
    operator_sigma_y(y2, qs, 1);
    y = qdpack_operator_alloc(qs);
    ws1 = qdpack_operator_alloc(qs);
    ws2 = qdpack_operator_alloc(qs);
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), y1, y2, (gsl_complex_rect(0.0,0.0)), y);
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), rho, y, (gsl_complex_rect(0.0,0.0)), ws1);
    qdpack_operator_blas_zgemm(CblasConjTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), rho, y, (gsl_complex_rect(0.0,0.0)), ws2);
    qdpack_operator_blas_zgemm(CblasNoTrans, CblasNoTrans, (gsl_complex_rect(1.0,0.0)), ws1, ws2, (gsl_complex_rect(0.0,0.0)), y);
    printf("rr = \n"); qdpack_operator_print(y);
    qdpack_matrix_eigen_zgeev(y->data, eval, evec, 0);
    {
        int i; int j; int n = 4;
        printf("Eigenvalues = \n");
        for (i = 0; i < n; i++)
        {
                struct gsl_complex z;
                z = qdpack_matrix_get(eval, i, 0);
                printf("\t(%f, %f)", ((z).dat[0]), ((z).dat[1]));
        }
        printf("\n");
        printf("evec = \n");
        for (i = 0; i < n; i++)
        {
                for (j = 0; j < n; j++)
                {
                        struct gsl_complex z;
                        z = qdpack_matrix_get(evec, i, j);
                        printf("\t(%f, %f)", ((z).dat[0]), ((z).dat[1]));
                }
                printf("\n");
        }
    }
    printf("c: eigvals: (%f, %f), (%f, %f), (%f, %f), (%f, %f)\n",
           ((qdpack_matrix_get(eval, 0, 0)).dat[0]), ((qdpack_matrix_get(eval, 0, 0)).dat[1]),
           ((qdpack_matrix_get(eval, 1, 0)).dat[0]), ((qdpack_matrix_get(eval, 1, 0)).dat[1]),
           ((qdpack_matrix_get(eval, 2, 0)).dat[0]), ((qdpack_matrix_get(eval, 2, 0)).dat[1]),
           ((qdpack_matrix_get(eval, 3, 0)).dat[0]), ((qdpack_matrix_get(eval, 3, 0)).dat[1]));
    c = + ((qdpack_matrix_get(eval, 0, 0)).dat[0]) - ((qdpack_matrix_get(eval, 1, 0)).dat[0])
        - ((qdpack_matrix_get(eval, 2, 0)).dat[0]) - ((qdpack_matrix_get(eval, 3, 0)).dat[0]);
    c = (c < 0.0) ? 0.0 : c;
    qdpack_operator_free(y1);
    qdpack_operator_free(y2);
    qdpack_operator_free(y);
    qdpack_operator_free(ws1);
    qdpack_operator_free(ws2);
    qdpack_matrix_free(evec);
    qdpack_matrix_free(eval);
    return c;
}

variable_types: struct qdpack_hilbert_space_t { int nsubsys; int *nstates; } *

variable_types: struct _qdpack_operator_t { int n; struct qdpack_hilbert_space_t *qs; struct _qdpack_matrix { int n; struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *data; } *data; } *

variable_types: c double

variable_types: j int

variable_types: n int

variable_types:
z
struct gsl_complex { double dat[2]; }


decompiled_code:
__int64 __fastcall gsl_ext_matrix_convert(_QWORD *a1, __int64 *a2)
{
__int64 result; // rax
unsigned int v3; // er15
__int64 v4; // r12
unsigned int v5; // ebp
__int64 i; // rbx
double v7; // xmm0_8

result = *a2;
if ( *a2 != *a1 )
return 0xFFFFFFFFLL;
if ( a2[1] != a1[1] )
return 0xFFFFFFFFLL;
if ( result )
{
v3 = 0;
v4 = 0LL;
v5 = 0;
for ( i = 0LL; ; i = 0LL )
{
do
{
v7 = gsl_matrix_get(a2, v4, i);
gsl_matrix_complex_set(a1, v4, i, v7, 0.0);
i = ++v5;
}
while ( v5 < (unsigned __int64)*a2 );
v4 = ++v3;
if ( *a2 <= (unsigned __int64)v3 )
break;
v5 = 0;
}
result = 0LL;
}
return result;
}

canonical_decompiled_code:
__int64 __fastcall func0(_QWORD *a1, __int64 *a2)
{
__int64 result; // rax
unsigned int v3; // er15
__int64 v4; // r12
unsigned int v5; // ebp
__int64 i; // rbx
double v7; // xmm0_8

result = *a2;
if ( *a2 != *a1 )
return 0xFFFFFFFFLL;
if ( a2[1] != a1[1] )
return 0xFFFFFFFFLL;
if ( result )
{
v3 = 0;
v4 = 0LL;
v5 = 0;
for ( i = 0LL; ; i = 0LL )
{
do
{
v7 = gsl_matrix_get(a2, v4, i);
gsl_matrix_complex_set(a1, v4, i, v7, 0.0);
i = ++v5;
}
while ( v5 < (unsigned __int64)*a2 );
v4 = ++v3;
if ( *a2 <= (unsigned __int64)v3 )
break;
v5 = 0;
}
result = 0LL;
}
return result;
}

original_code:
int
gsl_ext_matrix_convert(gsl_matrix_complex *dst, gsl_matrix *src)
{
    unsigned int i, j;
    gsl_complex z;
    if (src->size1 != dst->size1 || src->size2 != dst->size2)
        return -1;
    for (i = 0; i < src->size1; i++)
    {
        for (j = 0; j < src->size1; j++)
        {
            do {(&z)->dat[0]=(gsl_matrix_get(src, i, j)); (&z)->dat[1]=(0.0);} while(0);
            gsl_matrix_complex_set(dst, i, j, z);
        }
    }
    return 0;
}

canonical_original_code:
int gsl_ext_matrix_convert(struct gsl_matrix_complex *dst, struct gsl_matrix *src)
{
    unsigned int i; unsigned int j;
    struct gsl_complex z;
    if (src->size1 != dst->size1 || src->size2 != dst->size2)
        return -1;
    for (i = 0; i < src->size1; i++)
    {
        for (j = 0; j < src->size1; j++)
        {
            do {(&z)->dat[0]=(gsl_matrix_get(src, i, j)); (&z)->dat[1]=(0.0);} while(0);
            gsl_matrix_complex_set(dst, i, j, z);
        }
    }
    return 0;
}

variable_types: struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *

variable_types: struct gsl_matrix { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_struct { unsigned long size; double *data; } *block; int owner; } *

variable_types:
i
unsigned int

variable_types:
j
unsigned int

user_defined_types: struct gsl_matrix { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_struct *block; int owner; }


decompiled_code:
__int64 __fastcall gsl_ext_expm_complex(_QWORD *a1, __int64 a2)
{
__int64 v2; // r14
__int64 v3; // r12
int v4; // ebx

v2 = gsl_matrix_complex_alloc(*a1, a1[1]);
v3 = gsl_matrix_complex_alloc(*a1, a1[1]);
gsl_matrix_complex_set_identity(v2);
gsl_matrix_complex_set_identity(a2);
v4 = 50;
do
{
gsl_matrix_complex_scale(a2, 1.0 / (double)v4, 0.0);
gsl_blas_zgemm(111LL, 111LL, a2, a1, v3, 1.0, 0.0, 0.0, 0.0);
gsl_matrix_complex_add(v3, v2);
gsl_matrix_complex_memcpy(a2, v3);
--v4;
}
while ( v4 );
return 0LL;
}

canonical_decompiled_code:
__int64 __fastcall func1(_QWORD *a1, __int64 a2)
{
__int64 v2; // r14
__int64 v3; // r12
int v4; // ebx

v2 = gsl_matrix_complex_alloc(*a1, a1[1]);
v3 = gsl_matrix_complex_alloc(*a1, a1[1]);
gsl_matrix_complex_set_identity(v2);
gsl_matrix_complex_set_identity(a2);
v4 = 50;
do
{
gsl_matrix_complex_scale(a2, 1.0 / (double)v4, 0.0);
gsl_blas_zgemm(111LL, 111LL, a2, a1, v3, 1.0, 0.0, 0.0, 0.0);
gsl_matrix_complex_add(v3, v2);
gsl_matrix_complex_memcpy(a2, v3);
--v4;
}
while ( v4 );
return 0LL;
}

original_code:
int
gsl_ext_expm_complex(gsl_matrix_complex *A, gsl_matrix_complex *eA)
{
    int i;
    gsl_complex alpha, beta, z;
    gsl_matrix_complex *I, *T;
    I = gsl_matrix_complex_alloc(A->size1, A->size2);
    T = gsl_matrix_complex_alloc(A->size1, A->size2);
    do {(&alpha)->dat[0]=(1.0); (&alpha)->dat[1]=(0.0);} while(0);
    do {(&beta)->dat[0]=(0.0); (&beta)->dat[1]=(0.0);} while(0);
    gsl_matrix_complex_set_identity(I);
    gsl_matrix_complex_set_identity(eA);
    for (i = 50; i > 0; i--)
    {
        do {(&z)->dat[0]=(1.0 / i); (&z)->dat[1]=(0.0);} while(0);
        gsl_matrix_complex_scale(eA, z);
        gsl_blas_zgemm(CblasNoTrans, CblasNoTrans, alpha, eA, A, beta, T);
        gsl_matrix_complex_add(T, I);
        gsl_matrix_complex_memcpy(eA, T);
    }
    return 0;
}

canonical_original_code:
int gsl_ext_expm_complex(struct gsl_matrix_complex *A, struct gsl_matrix_complex *eA)
{
    int i;
    struct gsl_complex alpha; struct gsl_complex beta; struct gsl_complex z;
    struct gsl_matrix_complex *I; struct gsl_matrix_complex *T;
    I = gsl_matrix_complex_alloc(A->size1, A->size2);
    T = gsl_matrix_complex_alloc(A->size1, A->size2);
    do {(&alpha)->dat[0]=(1.0); (&alpha)->dat[1]=(0.0);} while(0);
    do {(&beta)->dat[0]=(0.0); (&beta)->dat[1]=(0.0);} while(0);
    gsl_matrix_complex_set_identity(I);
    gsl_matrix_complex_set_identity(eA);
    for (i = 50; i > 0; i--)
    {
        do {(&z)->dat[0]=(1.0 / i); (&z)->dat[1]=(0.0);} while(0);
        gsl_matrix_complex_scale(eA, z);
        gsl_blas_zgemm(CblasNoTrans, CblasNoTrans, alpha, eA, A, beta, T);
        gsl_matrix_complex_add(T, I);
        gsl_matrix_complex_memcpy(eA, T);
    }
    return 0;
}

variable_types:
eA
struct gsl_matrix_complex *

variable_types:
alpha
struct gsl_complex { double dat[2]; }

variable_types:
beta
struct gsl_complex { double dat[2]; }


original_code:
int
gsl_ext_eigen_zgeev(gsl_matrix_complex *A_gsl, gsl_matrix_complex *evec, gsl_vector_complex *eval)
{
    integer n,i,j,info,lwork, ldvl, ldvr, lda;
    doublecomplex *A,*vr,*vl,*w,*work;
    doublereal *rwork;
    char jobvl,jobvr;
    n = A_gsl->size1;
    A = (doublecomplex *)malloc((size_t)n * n * sizeof(doublecomplex));
    w = (doublecomplex *)malloc((size_t)n * sizeof(doublecomplex));
    vr = (doublecomplex *)malloc((size_t)n * n * sizeof(doublecomplex));
    vl = (doublecomplex *)malloc((size_t)n * n * sizeof(doublecomplex));
    lwork = 16 * n;
    work = (doublecomplex *)malloc((size_t)lwork * sizeof(doublecomplex));
    rwork = (doublereal *)malloc((size_t)lwork * sizeof(doublereal));
    for (i = 0; i < n; i++)
    {
        for (j = 0; j < n; j++)
        {
            gsl_complex z;
            double re,im;
            z = gsl_matrix_complex_get(A_gsl, i, j);
            re = ((z).dat[0]);
            im = ((z).dat[1]);
            A[j*n+i] = (doublecomplex){re,im};
        }
    }
    jobvl='N';
    jobvr='V';
    lda = n;
    ldvr = n;
    ldvl = n;
    zgeev_(&jobvl, &jobvr, &n, A, &lda, w, vl, &ldvl, vr, &ldvr, work, &lwork, rwork, &info);
    for (i = 0; i < n; i++)
    {
        gsl_complex z;
        do {(&z)->dat[0]=(w[i].r); (&z)->dat[1]=(w[i].i);} while(0);
        gsl_vector_complex_set(eval, i, z);
    }
    for (j = 0; j < n; j++)
    {
        for (i = 0; i < n; i++)
        {
            gsl_complex z;
            do {(&z)->dat[0]=(vr[j*n+i].r); (&z)->dat[1]=(vr[j*n+i].i);} while(0);
            gsl_matrix_complex_set(evec, i, j, z);
        }
    }
    if (info != 0)
    {
        printf("zgeev_: error: info = %d\n", (int)info);
    }
    free(A);
    free(w);
    free(vr);
    free(vl);
    free(work);
    free(rwork);
    return 0;
}

canonical_original_code:
int gsl_ext_eigen_zgeev(struct gsl_matrix_complex *A_gsl, struct gsl_matrix_complex *evec, struct gsl_vector_complex *eval)
{
    int n; int i; int j; int info; int lwork; int ldvl; int ldvr; int lda;
    struct doublecomplex *A; struct doublecomplex *vr; struct doublecomplex *vl; struct doublecomplex *w; struct doublecomplex *work;
    double *rwork;
    char jobvl; char jobvr;
    n = A_gsl->size1;
    A = (struct doublecomplex *)malloc((unsigned long)n * n * sizeof(doublecomplex));
    w = (struct doublecomplex *)malloc((unsigned long)n * sizeof(doublecomplex));
    vr = (struct doublecomplex *)malloc((unsigned long)n * n * sizeof(doublecomplex));
    vl = (struct doublecomplex *)malloc((unsigned long)n * n * sizeof(doublecomplex));
    lwork = 16 * n;
    work = (struct doublecomplex *)malloc((unsigned long)lwork * sizeof(doublecomplex));
    rwork = (double *)malloc((unsigned long)lwork * sizeof(doublereal));
    for (i = 0; i < n; i++)
    {
        for (j = 0; j < n; j++)
        {
            struct gsl_complex z;
            double re; double im;
            z = gsl_matrix_complex_get(A_gsl, i, j);
            re = ((z).dat[0]);
            im = ((z).dat[1]);
            A[j*n+i] = (struct doublecomplex){re,im};
        }
    }
    jobvl='N';
    jobvr='V';
    lda = n;
    ldvr = n;
    ldvl = n;
    zgeev_(&jobvl, &jobvr, &n, A, &lda, w, vl, &ldvl, vr, &ldvr, work, &lwork, rwork, &info);
    for (i = 0; i < n; i++)
    {
        struct gsl_complex z;
        do {(&z)->dat[0]=(w[i].r); (&z)->dat[1]=(w[i].i);} while(0);
        gsl_vector_complex_set(eval, i, z);
    }
    for (j = 0; j < n; j++)
    {
        for (i = 0; i < n; i++)
        {
            struct gsl_complex z;
            do {(&z)->dat[0]=(vr[j*n+i].r); (&z)->dat[1]=(vr[j*n+i].i);} while(0);
            gsl_matrix_complex_set(evec, i, j, z);
        }
    }
    if (info != 0)
    {
        printf("zgeev_: error: info = %d\n", (int)info);
    }
    free(A);
    free(w);
    free(vr);
    free(vl);
    free(work);
    free(rwork);
    return 0;
}

variable_types:
evec
struct gsl_matrix_complex *

variable_types: struct gsl_vector_complex { unsigned long size; unsigned long stride; double *data; struct gsl_block_complex_struct *block; int owner; } *

variable_types: i int

variable_types: info int

variable_types: lwork int

variable_types: ldvl int

variable_types: ldvr int

variable_types: lda int

variable_types:
A
struct doublecomplex { double r; } *

variable_types:
vr
struct doublecomplex { double r; } *

variable_types:
vl
struct doublecomplex { double r; } *

variable_types:
w
struct doublecomplex { double r; } *

variable_types:
work
struct doublecomplex { double r; } *

variable_types:
rwork
double *

variable_types: jobvl char

variable_types: jobvr char

variable_types: re double

variable_types: im double


original_code:
int
gsl_ext_eigen_sort(gsl_matrix_complex *evec, gsl_vector_complex *eval, int sort_order)
{
    gsl_complex z;
      gsl_matrix_complex *evec_copy;
    gsl_vector_complex *eval_copy;
    int *idx_map, i, j, idx_temp;
    double p1, p2;
    if ((evec->size1 != evec->size2) || (evec->size1 != eval->size))
    {
        return -1;
    }
    evec_copy = gsl_matrix_complex_alloc(evec->size1, evec->size2);
    eval_copy = gsl_vector_complex_alloc(eval->size);
    idx_map = (int *)malloc(sizeof(int) * eval->size);
    gsl_matrix_complex_memcpy(evec_copy, evec);
    gsl_vector_complex_memcpy(eval_copy, eval);
    for (i = 0; i < eval->size; i++)
    {
        idx_map[i] = i;
    }
    for (i = 0; i < eval->size - 1; i++)
    {
        for (j = i+1; j < eval->size; j++)
        {
            idx_temp = -1;
            if (sort_order == 0)
            {
                p1 = gsl_complex_abs(gsl_vector_complex_get(eval, idx_map[i]));
                p2 = gsl_complex_abs(gsl_vector_complex_get(eval, idx_map[j]));
                if (p1 > p2)
                {
                    idx_temp = idx_map[i];
                }
            }
            if (sort_order == 3)
            {
                p1 = gsl_complex_arg(gsl_vector_complex_get(eval, idx_map[i]));
                p2 = gsl_complex_arg(gsl_vector_complex_get(eval, idx_map[j]));
                if (p1 > 3.14159265358979323846) p1 -= 2*3.14159265358979323846;
                if (p1 <= -3.14159265358979323846) p1 += 2*3.14159265358979323846;
                if (p2 > 3.14159265358979323846) p2 -= 2*3.14159265358979323846;
                if (p2 <= -3.14159265358979323846) p2 += 2*3.14159265358979323846;
                if (p2 < p1)
                {
                    idx_temp = idx_map[i];
                }
            }
            if (idx_temp != -1)
            {
                idx_map[i] = idx_map[j];
                idx_map[j] = idx_temp;
            }
        }
    }
    for (i = 0; i < eval->size; i++)
    {
        for (j = 0; j < eval->size; j++)
        {
            z = gsl_matrix_complex_get(evec_copy, idx_map[i], j);
            gsl_matrix_complex_set(evec, i, j, z);
        }
        z = gsl_vector_complex_get(eval_copy, idx_map[i]);
        gsl_vector_complex_set(eval, i, z);
    }
    gsl_matrix_complex_free(evec_copy);
    gsl_vector_complex_free(eval_copy);
    free(idx_map);
    return 0;
}

canonical_original_code:
int gsl_ext_eigen_sort(struct gsl_matrix_complex *evec, struct gsl_vector_complex *eval, int sort_order)
{
    struct gsl_complex z;
      struct gsl_matrix_complex *evec_copy;
    struct gsl_vector_complex *eval_copy;
    int *idx_map; int i; int j; int idx_temp;
    double p1; double p2;
    if ((evec->size1 != evec->size2) || (evec->size1 != eval->size))
    {
        return -1;
    }
    evec_copy = gsl_matrix_complex_alloc(evec->size1, evec->size2);
    eval_copy = gsl_vector_complex_alloc(eval->size);
    idx_map = (int *)malloc(sizeof(int) * eval->size);
    gsl_matrix_complex_memcpy(evec_copy, evec);
    gsl_vector_complex_memcpy(eval_copy, eval);
    for (i = 0; i < eval->size; i++)
    {
        idx_map[i] = i;
    }
    for (i = 0; i < eval->size - 1; i++)
    {
        for (j = i+1; j < eval->size; j++)
        {
            idx_temp = -1;
            if (sort_order == 0)
            {
                p1 = gsl_complex_abs(gsl_vector_complex_get(eval, idx_map[i]));
                p2 = gsl_complex_abs(gsl_vector_complex_get(eval, idx_map[j]));
                if (p1 > p2)
                {
                    idx_temp = idx_map[i];
                }
            }
            if (sort_order == 3)
            {
                p1 = gsl_complex_arg(gsl_vector_complex_get(eval, idx_map[i]));
                p2 = gsl_complex_arg(gsl_vector_complex_get(eval, idx_map[j]));
                if (p1 > 3.14159265358979323846) p1 -= 2*3.14159265358979323846;
                if (p1 <= -3.14159265358979323846) p1 += 2*3.14159265358979323846;
                if (p2 > 3.14159265358979323846) p2 -= 2*3.14159265358979323846;
                if (p2 <= -3.14159265358979323846) p2 += 2*3.14159265358979323846;
                if (p2 < p1)
                {
                    idx_temp = idx_map[i];
                }
            }
            if (idx_temp != -1)
            {
                idx_map[i] = idx_map[j];
                idx_map[j] = idx_temp;
            }
        }
    }
    for (i = 0; i < eval->size; i++)
    {
        for (j = 0; j < eval->size; j++)
        {
            z = gsl_matrix_complex_get(evec_copy, idx_map[i], j);
            gsl_matrix_complex_set(evec, i, j, z);
        }
        z = gsl_vector_complex_get(eval_copy, idx_map[i]);
        gsl_vector_complex_set(eval, i, z);
    }
    gsl_matrix_complex_free(evec_copy);
    gsl_vector_complex_free(eval_copy);
    free(idx_map);
    return 0;
}

variable_types: struct gsl_matrix_complex { unsigned long size1; unsigned long size2; unsigned long tda; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *

variable_types: struct gsl_vector_complex { unsigned long size; unsigned long stride; double *data; struct gsl_block_complex_struct *block; int owner; } *

variable_types:
z
struct gsl_complex { double dat[2]; }

variable_types: struct gsl_vector_complex { unsigned long size; unsigned long stride; double *data; struct gsl_block_complex_struct { unsigned long size; double *data; } *block; int owner; } *

variable_types: j int