Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

#2331 fix bug in size of yS numpy array when simulation stops early #2337

Merged
merged 6 commits into from
Oct 6, 2022
Merged

#2331 fix bug in size of yS numpy array when simulation stops early #2337

Changes from 1 commit
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
128e02e
#2331 fix bug in size of yS numpy array when simulation stops early
martinjrobins Oct 4, 2022
7f08623
#2331 add test for sens with events, fix bug in python-based idaklu s…
martinjrobins Oct 6, 2022
fe0fb8b
style: pre-commit fixes
pre-commit-ci[bot] Oct 6, 2022
a6bb115
#2331 add changlog entry
martinjrobins Oct 6, 2022
25d1d0b
style: pre-commit fixes
pre-commit-ci[bot] Oct 6, 2022
749019d
Update CHANGELOG.md
valentinsulzer Oct 6, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
21 changes: 14 additions & 7 deletions pybamm/solvers/c_solvers/idaklu_casadi.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -373,6 +373,14 @@ int sensitivities_casadi(int Ns, realtype t, N_Vector yy, N_Vector yp,
// resvalsS now has (∂F/∂p i )
p_python_functions->sens();

//for (int i = 0; i < np; i++) {
// std::cout << "dF/dp[" << i << "] = [";
// for (int j = 0; j < p_python_functions->number_of_states; j++) {
// std::cout << NV_DATA_S(resvalS[i])[j] << " ";
// }
// std::cout << "]" << std::endl;
//}

martinjrobins marked this conversation as resolved.
Show resolved Hide resolved
for (int i = 0; i < np; i++) {
//std::cout << "dF/dp = [" << i << "] = [";
//for (int j = 0; j < p_python_functions->number_of_states; j++) {
Expand Down Expand Up @@ -457,7 +465,8 @@ Solution solve_casadi(np_array t_np, np_array y0_np, np_array yp0_np,
void *ida_mem; // pointer to memory
N_Vector yy, yp, avtol; // y, y', and absolute tolerance
N_Vector *yyS, *ypS; // y, y' for sensitivities
realtype rtol, *yval, *ypval, *atval, *ySval;
realtype rtol, *yval, *ypval, *atval;
std::vector<realtype *> ySval(number_of_parameters);
int retval;
SUNMatrix J;
SUNLinearSolver LS;
Expand All @@ -474,9 +483,6 @@ Solution solve_casadi(np_array t_np, np_array y0_np, np_array yp0_np,

// set initial value
yval = N_VGetArrayPointer(yy);
if (number_of_parameters > 0) {
ySval = N_VGetArrayPointer(yyS[0]);
}
ypval = N_VGetArrayPointer(yp);
atval = N_VGetArrayPointer(avtol);
int i;
Expand All @@ -488,6 +494,7 @@ Solution solve_casadi(np_array t_np, np_array y0_np, np_array yp0_np,
}

for (int is = 0 ; is < number_of_parameters; is++) {
ySval[is] = N_VGetArrayPointer(yyS[is]);
N_VConst(RCONST(0.0), yyS[is]);
N_VConst(RCONST(0.0), ypS[is]);
}
Expand Down Expand Up @@ -609,7 +616,7 @@ Solution solve_casadi(np_array t_np, np_array y0_np, np_array yp0_np,
for (int j = 0; j < number_of_parameters; j++) {
const int base_index = j * number_of_timesteps * number_of_states;
for (int k = 0; k < number_of_states; k++) {
yS_return[base_index + k] = ySval[j * number_of_states + k];
yS_return[base_index + k] = ySval[j][k];
}
}

Expand Down Expand Up @@ -649,7 +656,7 @@ Solution solve_casadi(np_array t_np, np_array y0_np, np_array yp0_np,
const int base_index = j * number_of_timesteps * number_of_states
+ t_i * number_of_states;
for (int k = 0; k < number_of_states; k++) {
yS_return[base_index + k] = ySval[j * number_of_states + k];
yS_return[base_index + k] = ySval[j][k];
}
}
t_i += 1;
Expand All @@ -665,7 +672,7 @@ Solution solve_casadi(np_array t_np, np_array y0_np, np_array yp0_np,
np_array t_ret = np_array(t_i, &t_return[0], free_t_when_done);
np_array y_ret = np_array(t_i * number_of_states, &y_return[0], free_y_when_done);
np_array yS_ret = np_array(
std::vector<ptrdiff_t>{number_of_parameters, t_i, number_of_states},
std::vector<ptrdiff_t>{number_of_parameters, number_of_timesteps, number_of_states},
&yS_return[0], free_yS_when_done
);

Expand Down