Write a function to implement Heun's method with an option to iterate the corrector equation to a specified stopping criterion. Your function need only work for a single ODE $($not a system$)$ and should accept the following inputs $($in order$)$:

A function dydt that defines the ODE to be solved. Your Heun function should accommodate optional parameter inputs to the dydt function.

An evenly spaced time vector that defines the step size, h for the integration $($by its increment$)$ and the time span $($first and last values$)$ over which the IVP is to be integrated.

The initial condition for y$.$

A stopping criterion for the corrector iteration. If the stopping criterion is left blank, your function should default to implementing Heun's method without iterating the corrector equation.

Varargin to accept any optional parameters needed to evaluate the dydt function.

Your function should output the following $($in order$)$:

A column vector of y values corresponding to the input time vector.

Note the test suite will test your function for the following two ODEs:

ODEs to be tested

This second ODE models constrained population growth. The parameter k$\_$gm is the maximum growth rate $($under unconstrained conditions$)$ and the parameter p$\_$max is the carrying capacity $($i$.$e$.$ maximum population$).$ The initial conditions, parameter values, and time span used in the test suite model growth of the Earth's human population from $1950-2100.$