Changeset 9082


Ignore:
Timestamp:
Jul 18, 2016 10:21:49 AM (3 years ago)
Author:
Christian Andersson
Message:

Cleaned up computation of the scaled vector norm. Note though that this should be moved to jmi_maths in the future. Related to ticket:4815

Location:
trunk/RuntimeLibrary/src/jmi
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/RuntimeLibrary/src/jmi/jmi_block_solver.c

    r9070 r9082  
    977977}
    978978
    979 double calculate_scaled_residual_norm(double* residual, double *f_scale, int n) {
    980     double scaledRes, norm = 0;
    981     int i;
    982     for(i=0; i<n; i++) {
    983         scaledRes = residual[i]*f_scale[i];
    984         norm = norm>= RAbs(scaledRes)?norm:RAbs(scaledRes);
    985     }
    986     return norm;
    987 }
    988 
     979int jmi_scaled_vector_norm(jmi_real_t *x, jmi_real_t *scale, jmi_int_t n, jmi_int_t NORM, jmi_real_t* out) {
     980    int ef = 0;
     981   
     982    if (NORM == JMI_NORM_MAX) {
     983        int i;
     984        jmi_real_t tmp;
     985        *out = 0.0;
     986        for (i = 0; i < n; i++) {
     987            tmp = JMI_ABS(x[i]*scale[i]);
     988            *out = *out >= tmp ? *out : tmp;
     989        }
     990    } else {
     991        ef = -1;
     992    }
     993   
     994    return ef;
     995}
     996
  • trunk/RuntimeLibrary/src/jmi/jmi_block_solver.h

    r9043 r9082  
    6565#define JMI_VAR_NOT_USED(x) ((void)x)
    6666
     67#define JMI_NORM_MAX 0
     68
    6769/** \brief Jacobian variability for the linear solver */
    6870typedef enum jmi_block_solver_jac_variability_t {
     
    395397void jmi_setup_f_residual_scaling(jmi_block_solver_t *block);
    396398
    397 /** \brief Calculate scaled residual max norm */
    398 double calculate_scaled_residual_norm(double* residual, double *f_scale, int n);
     399/** \brief Computes the scaled vector norm */
     400int jmi_scaled_vector_norm(jmi_real_t *x, jmi_real_t *scale, jmi_int_t n, jmi_int_t NORM, jmi_real_t* out);
    399401
    400402/** \brief Check and log illegal iv inputs */
  • trunk/RuntimeLibrary/src/jmi/jmi_linear_solver.c

    r9070 r9082  
    483483    /* Check if the calculated solution from minimum norm is a valid solution to the original problem */
    484484    if(solver->singular_jacobian==1) {
    485         double scaledMaxNorm;
     485        jmi_real_t scaled_max_norm;
     486        int ef = 0;
    486487        jmi_update_f_scale(block);
    487         scaledMaxNorm = calculate_scaled_residual_norm(solver->rhs, N_VGetArrayPointer(block->f_scale), block->n);
    488         if(scaledMaxNorm <= block->options->res_tol) {
     488       
     489        ef = jmi_scaled_vector_norm(solver->rhs, N_VGetArrayPointer(block->f_scale), block->n, JMI_NORM_MAX, &scaled_max_norm);
     490        if (ef == -1) {
     491            jmi_log_node(block->log, logError, "NormFailure", "Failed to compute the scaled residual norm to the linear system in <block: %s>", block->label);
     492        }
     493       
     494        if(scaled_max_norm <= block->options->res_tol) {
    489495            if(block->callbacks->log_options.log_level >= 5){
    490496                jmi_log_node(block->log, logInfo, "Info", "Successfully calculated the minimum norm solution to the linear system in <block: %s>", block->label);
     
    494500            destnode = jmi_log_enter_fmt(block->log, logError, "UnsolveableLinearSystem", "Failed to calculate a valid minimum norm solution to the linear system in <block: %s> at <t: %f>", block->label, block->cur_time);
    495501            jmi_log_reals(block->log, destnode, logError, "residuals", solver->rhs, block->n);
    496             jmi_log_reals(block->log, destnode, logError, "scaled_max_norm", &(scaledMaxNorm), 1);
    497             jmi_log_reals(block->log, destnode, logError, "tolerance", &(block->options->res_tol), 1);
     502            jmi_log_reals(block->log, destnode, logError, "scaled_max_norm", &(scaled_max_norm), 1);
     503            jmi_log_reals(block->log, destnode, logError, "tolerance", &(block->options->res_tol), 1);
    498504            jmi_log_leave(block->log, destnode);
    499505        }
Note: See TracChangeset for help on using the changeset viewer.