% Applies indexed equations to all nodes in a domain 
function dT = IndexFDEqs(Domain,T,t)
  if iscell(T)             % multi domain problem T=T{j}(m,n)
    if size(Domain{1},2)==1
      error('Domain must use cell format Domain{j}');
    end
    for j=1:length(T)
      dT{j}=LoopEqs(Domain{j},T,t,size(T{j}),j);
    end
  else                     % single domain problem T=T(m,n)
    dT=LoopEqs(Domain,T,t,size(T),0);
  end
end

function dT = LoopEqs(Domain,T,t,DimT,j)
  [FDEqs,nodes]=Domain{1:2};
  dT=zeros(DimT);
  Dim=sum(DimT>1);
  for i=1:length(FDEqs)
    if Dim<=1                             % 1-D
      Nods={nodes{i}};
    elseif Dim==2                         % 2-D
      Nods={nodes{i,1},nodes{i,2}};
    elseif Dim==3                         % 3-D
      Nods={nodes{i,1},nodes{i,2},nodes{i,3}};
    end
    if j
      FunIndx={Nods{:},j};
    else
      FunIndx=Nods;
    end
    dT(Nods{:})=FDEqs{i}(T,t,FunIndx{:});
  end
end
