diff --git a/ProcessLib/PhaseFieldAcid/PhaseFieldAcidFEM.h b/ProcessLib/PhaseFieldAcid/PhaseFieldAcidFEM.h index 4a08c2ff23945a8218fc9fef670c453b7488c59e..1cf19da9b2ceeb6c912ab615106b0bad76a2393f 100644 --- a/ProcessLib/PhaseFieldAcid/PhaseFieldAcidFEM.h +++ b/ProcessLib/PhaseFieldAcid/PhaseFieldAcidFEM.h @@ -277,16 +277,17 @@ private: double const alpha = _process_data.alpha(t, x_position)[0]; double const rrc = _process_data.rrc(t, x_position)[0]; double const da = rrc * epsilon / D; - double const tau = 1.0;//epsilon * epsilon / D; + double const tau = 1.0; // epsilon * epsilon / D; // lambda = // tau * D / epsilon / epsilon / (alpha * (5.0 / 3.0 + sqrt(2.) / da)) - double const lambda = 1.2;// / (alpha * (5.0 / 3.0 + sqrt(2.) / da)); + double const lambda = 1.2; // / (alpha * (5.0 / 3.0 + sqrt(2.) / da)); double const phi_avg = ph.sum() / ph.size(); double const c_avg = c.sum() / c.size(); (*_process_data.element_f)[_element.getID()] = f(phi_avg, c_avg, lambda); GlobalDimVectorType v_avg = GlobalDimVectorType::Zero(GlobalDim); + double element_kappa = (*_process_data.element_kappa)[_element.getID()]; for (int ip = 0; ip < n_integration_points; ip++) { @@ -330,8 +331,10 @@ private: // local_b.noalias() += // (N * psi_ip.dot(v_ip) + v_ip.transpose() * dNdx) * epsilon * // epsilon * w; + local_b.noalias() -= - kappa_ip * epsilon * epsilon * v_ip.norm() * N * w; + kappa_ip * epsilon * epsilon * v_ip.norm() * N * w; + // element_kappa * epsilon * epsilon * v_avg.norm()*N*w; if (_element.getID() == -1) // && ip == 0) { @@ -422,12 +425,13 @@ private: GlobalDimMatrixType const grad_v_ip = dNdx * v_at_nodes.transpose(); double const div_grad_phi_ip = grad_v_ip.diagonal().sum(); - double source = 0.0; + double source = 1.0; // TODO (naumov) Add prj file switch to disable it. - if (grad_ph_norm > std::numeric_limits<double>::epsilon()) - { - source = 0.0*(D * div_grad_phi_ip - ph_dot_ip) / rrc / grad_ph_norm; - } +// if (grad_ph_norm > std::numeric_limits<double>::epsilon()) +// { +// source += 0.0 * (D * div_grad_phi_ip - ph_dot_ip) / rrc / +// grad_ph_norm; +// } local_M.noalias() += w * N.transpose() * N;