     gsl_matrix *Q = gsl_matrix_alloc(p_basis->num, p_basis->num);
     form_Q(p_basis, Q);
     */
+
+    t1 = clock();
+    time_in_sec = (t1 - t0) / (double)CLOCKS_PER_SEC;
+    time_txt += "Time_SCF_Init = " + std::to_string(time_in_sec) + " sec\n";
+    time_total += time_in_sec;
 
 
     // start SCF iterations

            */
        }
    }


    // count number of primitive integrals in a <bra| or |ket>

    // SCF converged
    fprintf(stdout, "SCF converged! E_total = %20.10f\n", ene_total);

    t0 = clock();
    time_in_sec = (t0 - t1) / (double)CLOCKS_PER_SEC;
    time_txt += "Time_SCF_Conv = " + std::to_string(time_in_sec) + " sec\n";
    time_total += time_in_sec;


    // print MO information
    fprintf(stdout, "%5s %10s %15s %12s\n", "MO", "State", "E(Eh)", "E(eV)");
    for (ibasis = 0; ibasis < p_basis->num; ++ ibasis)
