function [ctempnew] = matsolveD(ML,MR,ctempold,bcvector);
%%
%% Solve ML*ctempnew = MR*ctempold + (lambda/2) c(x=0,old) + (lambda/2) c(x=0,new).
%% bcvector holds (lambda/2) c(x=0,old) + (lambda/2) c(x=0,new).
RHS = MR*ctempold + bcvector;
ctempnew = ML\RHS;
return;