Doing this with anonymous functions is not the way to go (not to mention, your last two functions, dc2D_dt and dcD_dt, have t as an input, but t doesn't actually show up in the function anywhere).
First, write down (i.e., by hand, in regular math notation) your set of equations in the form f_i(x) = 0, where x is your vector of unknowns, and i = 1,2,...,7 indexes your 7 different equations.
Next, you'll want to write a regular MATLAB function (i.e., not an anonymous function) that takes x as the input and returns the vector f(x) = [f_1(x); f_2(x); ..., f_7(x)] as its output.
Then you'll need to choose a MATLAB solver. There are various options available that each have pros and cons. Programatically, you need to choose between a problem-based approach and a solver-based one. From there, you then need to choose a particular solver.
The most widely applicable of the options is probably the solver-based fsolve (though depending on the structure of your problem, and your preferred programming approach, there may be something else that's preferable). fsolve will attemp to find an x that makes each of the elements of f(x) simultaneously as close to zero as possible.
