Though my title is quite general (please feel free to edit), indeed, I wonder if the following models can be solved in Gurobi or CPLEX. Model 2 is just an alternative one to Model 1. Although we discussed about an MIP reformulation applicability in here for a similar version of the problem, I would like to directly solve the problem as is. If Gurobi or CPLEX cannot handle it, can you please suggest some other solvers that can handle 10,000 variables as $|P|\approx10000$. Note that all notations other than $h_p$ (which are continuous in $\mathbb{R}^+$) represent some non-negative constant.
Model 1:
\begin{alignat}2\max &\quad \sum_{\substack{p\in P,\\s\in S_p}}e_{ps}\left(H^+-h_p\right)\tag1\\\text{s.t.}&\quad \sum_{p\in P}\left(\frac{k_p\tau_p}{h_p}-\sum_{s\in S_p}e_{ps}\right)\leq \kappa\tag2\\&\quad\sum_{\{p\in P|d_p=t,\theta_p=j\}}\frac{\tau_p}{h_p} \leq M_{tj} \qquad \forall t\in T, j\in J\tag3\\&\quad f_{ps}h_p \leq B_p \qquad \forall p\in P, s\in S_p\tag4\\&\quad h_p\in \mathbb{R}^+, H^- \leq h_p \leq H^+.\end{alignat}
Model 2:
\begin{alignat}2\max &\quad \sum_{\substack{p\in P,\\s\in S_p}}\frac{e_{ps}}{2}\left(30-\frac{30}{h_p}\right)\tag1\\\text{s.t.}&\quad \sum_{p\in P}\left(\frac{k_p\tau_ph_p}{60}-\sum_{s\in S_p}e_{ps}\right)\leq \kappa\tag2\\&\quad\sum_{\{p\in P|d_p=t,\theta_p=j\}}\frac{\tau_ph_p}{30} \leq M_{tj} \qquad \forall t\in T, j\in J\tag3\\&\quad h_p\geq \frac{e_{ps}}{B_p} \qquad \forall p\in P, s\in S_p\tag4\\&\quad h_p\in \mathbb{R}^+, H^- \leq h_p \leq H^+.\end{alignat}
