diff --git a/sPuReMD/src/charges.c b/sPuReMD/src/charges.c
index 4cbaa806d25176b45d611c275a424b800472cbac..7561df494dd1246dd0e51cd8a812d2bd9e90be4a 100644
--- a/sPuReMD/src/charges.c
+++ b/sPuReMD/src/charges.c
@@ -314,18 +314,19 @@ static void Compute_Preconditioner_QEq( const reax_system * const system,
// break;
//
// case ICHOLT_PC:
-// data->timing.cm_solver_pre_comp +=
-// ICHOLT( Hptr, workspace->droptol, workspace->L_EE, workspace->U_EE );
+// fprintf( stderr, "[ERROR] ICHOLT is not supported for indefinite, symmetric matrices of EE. Terminating...\n" );
+// exit( INVALID_INPUT );
// break;
//
-// case ILUT_PC:
-// data->timing.cm_solver_pre_comp +=
-// ILUT( Hptr, workspace->droptol, workspace->L, workspace->U );
-// break;
+// case ILUT_PC:
+// data->timing.cm_solver_pre_comp +=
+// ILUT( Hptr, workspace->droptol, workspace->L, workspace->U );
+// break;
//
-// case ILUTP_PC:
-// data->timing.cm_solver_pre_comp +=
-// ILUTP( Hptr, workspace->droptol, workspace->L, workspace->U );
+// case ILUTP_PC:
+// data->timing.cm_solver_pre_comp +=
+// ILUTP( Hptr, workspace->droptol, workspace->L, workspace->U );
+// break;
//
// case FG_ILUT_PC:
// data->timing.cm_solver_pre_comp +=
@@ -538,8 +539,7 @@ static void Compute_Preconditioner_EE( const reax_system * const system,
Hptr = create_test_mat( );
#endif
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 1.0;
@@ -564,8 +564,8 @@ static void Compute_Preconditioner_EE( const reax_system * const system,
break;
case ICHOLT_PC:
- data->timing.cm_solver_pre_comp +=
- ICHOLT( Hptr, workspace->droptol, workspace->L, workspace->U );
+ fprintf( stderr, "[ERROR] ICHOLT is not supported for indefinite, symmetric matrices of EE. Terminating...\n" );
+ exit( INVALID_INPUT );
break;
case ILUT_PC:
@@ -576,6 +576,7 @@ static void Compute_Preconditioner_EE( const reax_system * const system,
case ILUTP_PC:
data->timing.cm_solver_pre_comp +=
ILUTP( Hptr, workspace->droptol, workspace->L, workspace->U );
+ break;
case FG_ILUT_PC:
data->timing.cm_solver_pre_comp +=
@@ -613,8 +614,7 @@ static void Compute_Preconditioner_EE( const reax_system * const system,
}
}
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 0.0;
@@ -663,8 +663,7 @@ static void Compute_Preconditioner_ACKS2( const reax_system * const system,
Hptr = create_test_mat( );
#endif
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 1.0;
@@ -690,8 +689,8 @@ static void Compute_Preconditioner_ACKS2( const reax_system * const system,
break;
case ICHOLT_PC:
- data->timing.cm_solver_pre_comp +=
- ICHOLT( Hptr, workspace->droptol, workspace->L, workspace->U );
+ fprintf( stderr, "[ERROR] ICHOLT is not supported for indefinite, symmetric matrices of ACKS2. Terminating...\n" );
+ exit( INVALID_INPUT );
break;
case ILUT_PC:
@@ -702,6 +701,7 @@ static void Compute_Preconditioner_ACKS2( const reax_system * const system,
case ILUTP_PC:
data->timing.cm_solver_pre_comp +=
ILUTP( Hptr, workspace->droptol, workspace->L, workspace->U );
+ break;
case FG_ILUT_PC:
data->timing.cm_solver_pre_comp +=
@@ -739,8 +739,7 @@ static void Compute_Preconditioner_ACKS2( const reax_system * const system,
}
}
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 0.0;
@@ -894,9 +893,7 @@ static void Setup_Preconditioner_EE( const reax_system * const system,
}
data->timing.cm_sort_mat_rows += Get_Timing_Info( time );
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
- || control->cm_solver_pre_comp_type == ILUTP_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 1.0;
@@ -916,23 +913,8 @@ static void Setup_Preconditioner_EE( const reax_system * const system,
break;
case ICHOLT_PC:
- Calculate_Droptol( Hptr, workspace->droptol, control->cm_solver_pre_comp_droptol );
-
- fillin = Estimate_LU_Fill( Hptr, workspace->droptol );
-
- if ( workspace->L == NULL )
- {
- Allocate_Matrix( &workspace->L, system->N_cm, fillin + system->N_cm );
- Allocate_Matrix( &workspace->U, system->N_cm, fillin + system->N_cm );
- }
- else if ( workspace->L->m < fillin + system->N_cm )
- {
- Deallocate_Matrix( workspace->L );
- Deallocate_Matrix( workspace->U );
-
- Allocate_Matrix( &workspace->L, system->N_cm, fillin + system->N_cm );
- Allocate_Matrix( &workspace->U, system->N_cm, fillin + system->N_cm );
- }
+ fprintf( stderr, "[ERROR] ICHOLT is not supported for indefinite, symmetric matrices of EE. Terminating...\n" );
+ exit( INVALID_INPUT );
break;
case ILUT_PC:
@@ -972,9 +954,7 @@ static void Setup_Preconditioner_EE( const reax_system * const system,
break;
}
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
- || control->cm_solver_pre_comp_type == ILUTP_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 0.0;
@@ -1010,9 +990,7 @@ static void Setup_Preconditioner_ACKS2( const reax_system * const system,
}
data->timing.cm_sort_mat_rows += Get_Timing_Info( time );
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
- || control->cm_solver_pre_comp_type == ILUTP_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 1.0;
@@ -1033,24 +1011,8 @@ static void Setup_Preconditioner_ACKS2( const reax_system * const system,
break;
case ICHOLT_PC:
- Calculate_Droptol( Hptr, workspace->droptol, control->cm_solver_pre_comp_droptol );
-
- fillin = Estimate_LU_Fill( Hptr, workspace->droptol );
-
- if ( workspace->L == NULL )
- {
- Allocate_Matrix( &workspace->L, Hptr->n, fillin );
- Allocate_Matrix( &workspace->U, Hptr->n, fillin );
- }
- else if ( workspace->L->m < fillin )
- {
- Deallocate_Matrix( workspace->L );
- Deallocate_Matrix( workspace->U );
-
- /* factors have sparsity pattern as H */
- Allocate_Matrix( &workspace->L, Hptr->n, fillin );
- Allocate_Matrix( &workspace->U, Hptr->n, fillin );
- }
+ fprintf( stderr, "[ERROR] ICHOLT is not supported for indefinite, symmetric matrices of ACKS2. Terminating...\n" );
+ exit( INVALID_INPUT );
break;
case ILUT_PC:
@@ -1089,9 +1051,7 @@ static void Setup_Preconditioner_ACKS2( const reax_system * const system,
break;
}
- if ( control->cm_solver_pre_comp_type == ICHOLT_PC
- || control->cm_solver_pre_comp_type == ILUT_PC
- || control->cm_solver_pre_comp_type == ILUTP_PC
+ if ( control->cm_solver_pre_comp_type == ILUT_PC
|| control->cm_solver_pre_comp_type == FG_ILUT_PC )
{
Hptr->val[Hptr->start[system->N + 1] - 1] = 0.0;
diff --git a/sPuReMD/src/geo_tools.c b/sPuReMD/src/geo_tools.c
index cbf9ca42b1b9b85f95fe0c0712f2d9989ddd1333..45e36abde14d6511fcca8fc467474ed0ff23b126 100644
--- a/sPuReMD/src/geo_tools.c
+++ b/sPuReMD/src/geo_tools.c
@@ -97,7 +97,7 @@ static int Read_Box_Info( reax_system *system, FILE *geo, int geo_format )
/* compute full volume tensor from the angles */
Setup_Box( atof(s_a), atof(s_b), atof(s_c),
atof(s_alpha), atof(s_beta), atof(s_gamma),
- &(system->box) );
+ &system->box );
ret = SUCCESS;
}
diff --git a/sPuReMD/src/print_utils.c b/sPuReMD/src/print_utils.c
index 302eef023e5890aa91785b8447afdcd2e006de88..e95ac23dc03aec6ccd75223fa1cdc3b1dfc05edc 100644
--- a/sPuReMD/src/print_utils.c
+++ b/sPuReMD/src/print_utils.c
@@ -867,6 +867,42 @@ void Print_Sparse_Matrix2( sparse_matrix *A, char *fname, char *mode )
}
+/* Read sparse matrix in COO format (1-based indexing)
+ * and store in a CSR sparse matrix (preallocated)
+ *
+ * Note: the file must be sorted in increasing order of
+ * columns than rows */
+void Read_Sparse_Matrix2( sparse_matrix *A, char *fname )
+{
+ int top, cur_row, row, col, val;
+ FILE *f;
+
+ f = sfopen( fname, "r" );
+ top = 0;
+ cur_row = 0;
+
+ A->start[cur_row] = top;
+
+ while ( fscanf( f, "%d %d %f", &row, &col, &val ) == 3 )
+ {
+ if ( cur_row != row - 1 )
+ {
+ cur_row++;
+ A->start[cur_row] = top;
+ }
+
+ A->j[top] = col - 1;
+ A->val[top] = val;
+
+ top++;
+ }
+
+ A->start[A->n] = top;
+
+ sfclose( f, "Read_Sparse_Matrix2::f" );
+}
+
+
/* Note: watch out for portability issues with endianness
* due to serialization of numeric types (integer, IEEE 754) */
void Print_Sparse_Matrix_Binary( sparse_matrix *A, char *fname )