function Start_Hier_Community(varargin)

multiple_population_handle = findobj('Tag','Fig_multiple_population');

param = get(multiple_population_handle,'UserData');

V = 15;

if(strcmp(param{1},'IDLE'))
   
   populations_handle = findobj(multiple_population_handle,'Tag','PopupMenu_populations');
   number_of_populations_handle = findobj(multiple_population_handle,'Tag','EditText_number_of_populations');
   maximum_time_handle = findobj(multiple_population_handle,'Tag','EditText_maximum_time');
   number_of_time_steps_handle = findobj(multiple_population_handle,'Tag','EditText_number_of_time_steps');
   maximum_x_handle = findobj(multiple_population_handle,'Tag','EditText_maximum_x');
   number_of_x_intervals_handle = findobj(multiple_population_handle,'Tag','EditText_number_of_x_intervals');
   status_bar_handle = findobj(multiple_population_handle,'Tag','StaticText_status_bar');
   start_handle = findobj(multiple_population_handle,'Tag','Pushbutton_start');
   axes_density_handle = findobj(multiple_population_handle,'Tag','Axes_population_density');
   axes_moment_handle = findobj(multiple_population_handle,'Tag','Axes_population_moment');
   cont_handle = findobj(multiple_population_handle,'Tag','Radiobutton_cont');
   order_handle = findobj(multiple_population_handle,'Tag','PopupMenu_order');
   
   set(status_bar_handle,'String','Initializing');
   
   num_pop = str2num(get(number_of_populations_handle,'String'));
   tmax = str2num(get(maximum_time_handle,'String'));
   tsteps = str2num(get(number_of_time_steps_handle,'String'));
   xmax = str2num(get(maximum_x_handle,'String'));
   xintv = str2num(get(number_of_x_intervals_handle,'String'));
   CONT = get(cont_handle,'Value');
   order_scheme = get(order_handle,'Value');
      
   param{1} = 'ACTIVE';
   set(multiple_population_handle,'UserData',param);
   pop_param = param{3};
   
   for i = 1:num_pop
      uo_f{i} = getfield(pop_param(i),'uo');
      g_f{i} = getfield(pop_param(i),'g');
      m_f{i} = getfield(pop_param(i),'m');
      b_f{i} = getfield(pop_param(i),'b');
      c_f{i} = getfield(pop_param(i),'c');
      w_f{i} = getfield(pop_param(i),'w');
      if(~CONT)
         pos = findstr(w_f{i},'<=');
         w_f{i}(pos+1)='';
      end
      Pvar{i} = ['P',num2str(i)];
   end
   
   for i = 1:num_pop
      [uo_var{i},uo_pos{i}] = Find_Var(uo_f{i},{'x'}); 
      [g_var{i},g_pos{i}] = Find_Var(g_f{i},[{'x'},Pvar]); 
      [m_var{i},m_pos{i}] = Find_Var(m_f{i},[{'x'},Pvar]); 
      [b_var{i},b_pos{i}] = Find_Var(b_f{i},[{'x'},Pvar]); 
      [c_var{i},c_pos{i}] = Find_Var(c_f{i},{'t'}); 
      [w_var{i},w_pos{i}] = Find_Var(w_f{i},{'x','y'}); 
   end
   
   dx = xmax/xintv;
   x=0:dx:xmax;
   for i=1:length(x)
       xs{i}=num2str(x(i),V);
   end
   
   dt = tmax/tsteps;
   t=0:dt:tmax;
   for i = 1:length(t)
       ts{i} = num2str(t(i),V);
   end
   t_end = tsteps;
   
   for i=1:num_pop
       hierarchy(i)=any(strcmp('x',w_var{i}));
   end
       
   if(isempty(varargin))
             
      t_begin = 1;
      for i = 1:num_pop
          
          for j = 1:xintv+1
              u{i}(1,j) = MyInline(uo_f{i},{'x'},xs(j),uo_var{i},uo_pos{i});
              
              if(hierarchy(i))
                 for n = 1:xintv+1
                     w{i}(j,n) = MyInline(w_f{i},{'x','y'},[xs(j),xs(n)],w_var{i},w_pos{i});
                 end
              else
                  w{i}(1:xintv+1,j)=MyInline(w_f{i},{'y'},xs(j),w_var{i},w_pos{i});
              end
          end
          
          for j = 1:xintv+1
              if (CONT)
                 tpop = dx*(sum(w{i}(j,2:end).*u{i}(1,2:end)));
              else
                 tpop = dx*(sum(w{i}(j,2:end).*u{i}(1,2:end))); 
              end
              Pin{i,j} = num2str(tpop,V);
              Ptotal{i}(t_begin,j)= tpop;
           end
       end
       
   else
      
      u = param{4}{3};
      Ptotal=param{4}{4};
      for i=1:num_pop
          for j=1:xintv+1
              if(hierarchy(i))
                  for n=1:xintv+1
                      w{i}(j,n)=MyInline(w_f{i},{'x','y'},[xs(j),xs(n)],w_var{i},w_pos{i});
                  end
              else
                  w{i}(1:xintv+1,j)=MyInline(w_f{i},{'y'},xs(j),w_var{i},w_pos{i});
              end
          end
          for j = 1:xintv+1
              if (CONT)
                 tpop = dx*(sum(w{i}(j,2:end).*u{i}(end,2:end)));
              else
                 tpop = dx*(sum(w{i}(j,2:end).*u{i}(end,2:end))); 
              end
              Pin{i,j} = num2str(tpop,V);
              Ptotal{i}(end,j)= tpop;
          end
      end
      t_begin = length(u{1}(:,1));
   end
   
   for k = t_begin : t_end
      
      set(status_bar_handle,'String',num2str(k));
      drawnow;
      
      param = get(multiple_population_handle,'UserData');
      if(strcmp(param{1},'IDLE'))
         param{1} = 'IDLE';
         current_population = get(populations_handle,'Value');
         num_stop=length(u{current_population}(:,1));
         param{4} = {x,t(1:num_stop),u,Ptotal};
         set(multiple_population_handle,'UserData',param);
   
         axes(axes_density_handle);
         cla;
         mesh(x,t(1:num_stop),u{current_population})
         title(['Population ',num2str(current_population),':  Population Density'],...
                  'Color','blue','FontSize',10,'Units','normalized');
         xlabel('x');
         ylabel('t');
         set(axes_density_handle,'Tag','Axes_population_density');
   
         axes(axes_moment_handle);
         cla;
         if (hierarchy(current_population))
             mesh(x,t(1:num_stop),Ptotal{current_population})
             xlabel('x');
             ylabel('t');
         else
             plot(t(1:num_stop),Ptotal{current_population}(:,1))
             xlabel('t');
         end
         title(['Population ',num2str(current_population),':  Total Population'],...
                  'Color','blue','FontSize',10,'Units','normalized');
         set(axes_moment_handle,'Tag','Axes_population_moment'); 
         
         set(status_bar_handle,'String','Simulation Stop By User');
         return;
      end
      
      for i = 1:num_pop
         
         c{i} = MyInline(c_f{i},{'t'},ts(k),c_var{i},c_pos{i});
         for j = 1:xintv+1
            g{i}(j) = MyInline(g_f{i},[{'x'},Pvar],[xs(j),Pin(:,j)'],g_var{i},g_pos{i});     
            b{i}(j) = MyInline(b_f{i},[{'x'},Pvar],[xs(j),Pin(:,j)'],b_var{i},b_pos{i});
            m{i}(j) = MyInline(m_f{i},[{'x'},Pvar],[xs(j),Pin(:,j)'],m_var{i},m_pos{i});
            d{i}(j) = 1 + dt*g{i}(j)/dx + dt*m{i}(j);
         end
         
         B(1,1)=g{i}(1);
         for j=1:xintv
            B(j+1,j)=-dt/dx*g{i}(j);
            B(j+1,j+1)=d{i}(j+1);
         end
         f=[c{i}+dx*sum(b{i}(2:end).*u{i}(k,2:end)),u{i}(k,2:end)]';
         u{i}(k+1,:)=(B\f)';
         
         for j = 1:xintv+1
            if(CONT)
               tpop = dx*(sum(w{i}(j,2:end).*u{i}(k+1,2:end)));
            else
               tpop = dx*(sum(w{i}(j,2:end).*u{i}(k+1,2:end))); 
            end
            Pin{i,j} = num2str(tpop,V);
            Ptotal{i}(k+1,j)= tpop;
         end
         
      end
   end
   
   if(order_scheme==2)
       x_half=0:dx/2:xmax;
       for i=1:length(x_half)
           xs{i}=num2str(x_half(i),V);
       end
       t_half=0:dt/2:tmax;
       for i = 1:length(t_half)
           ts{i} = num2str(t_half(i),V);
       end
       t_end = 2*tsteps;
   
       %if(isempty(varargin))
             
          t_begin = 1;
          for i = 1:num_pop
          
              for j = 1:2*xintv+1
                  u_half{i}(1,j) = MyInline(uo_f{i},{'x'},xs(j),uo_var{i},uo_pos{i});
              
                 if(hierarchy(i))
                    for n = 1:2*xintv+1
                        w{i}(j,n) = MyInline(w_f{i},{'x','y'},[xs(j),xs(n)],w_var{i},w_pos{i});
                    end
                 else
                    w{i}(1:2*xintv+1,j)=MyInline(w_f{i},{'y'},xs(j),w_var{i},w_pos{i});
                 end
              end
          
              for j = 1:2*xintv+1
                 if (CONT)
                    tpop = dx/2*(sum(w{i}(j,2:end).*u_half{i}(1,2:end)));
                 else
                    tpop = dx/2*(sum(w{i}(j,2:end).*u_half{i}(1,2:end))); 
                 end
                 Pin{i,j} = num2str(tpop,V);
                 Ptotal_half{i}(t_begin,j)= tpop;
             end
         end
         %end
      for k=t_begin:t_end
       
          set(status_bar_handle,'String',num2str(k));
          drawnow;
          
         for i = 1:num_pop
            c{i} = MyInline(c_f{i},{'t'},ts(k),c_var{i},c_pos{i});
            for j = 1:2*xintv+1
                g{i}(j) = MyInline(g_f{i},[{'x'},Pvar],[xs(j),Pin(:,j)'],g_var{i},g_pos{i});     
                b{i}(j) = MyInline(b_f{i},[{'x'},Pvar],[xs(j),Pin(:,j)'],b_var{i},b_pos{i});
                m{i}(j) = MyInline(m_f{i},[{'x'},Pvar],[xs(j),Pin(:,j)'],m_var{i},m_pos{i});
                d{i}(j) = 1 + dt*g{i}(j)/dx + dt/2*m{i}(j);
            end
         
            B(1,1)=g{i}(1);
            for j=1:2*xintv
                B(j+1,j)=-dt/dx*g{i}(j);
               B(j+1,j+1)=d{i}(j+1);
            end
            f=[c{i}+dx/2*sum(b{i}(2:end).*u_half{i}(k,2:end)),u_half{i}(k,2:end)]';
            u_half{i}(k+1,:)=(B\f)';
         
            for j = 1:2*xintv+1
               if(CONT)
                 tpop = dx/2*(sum(w{i}(j,2:end).*u_half{i}(k+1,2:end)));
               else
                 tpop = dx/2*(sum(w{i}(j,2:end).*u_half{i}(k+1,2:end))); 
               end
               Pin{i,j} = num2str(tpop,V);
               Ptotal_half{i}(k+1,j)= tpop;
            end
         
         end
      end 
      for i=1:num_pop
          for k=t_begin:t_end/2+1
              for j=1:xintv+1
                  u{i}(k,j)=1/3*(4*u_half{i}(2*k-1,2*j-1)-u{i}(k,j));
                  Ptotal{i}(k,j)=1/3*(4*Ptotal_half{i}(2*k-1,2*j-1)-Ptotal{i}(k,j));
              end
          end
      end
      
   end
   
   
   param{1} = 'IDLE';
   param{4} = {x,t,u,Ptotal};
   set(multiple_population_handle,'UserData',param);
   current_population = get(populations_handle,'Value');
   
   axes(axes_density_handle);
   cla;
   mesh(x,t,u{current_population})
   title(['Population ',num2str(current_population),':  Population Density'],...
      'Color','blue','FontSize',10,'Units','normalized');
   xlabel('x');
   ylabel('t');
   set(axes_density_handle,'Tag','Axes_population_density');
   
   axes(axes_moment_handle);
   cla;
   if (hierarchy(current_population))
       mesh(x,t,Ptotal{current_population})
       xlabel('x');
       ylabel('t');
   else
       plot(t,Ptotal{current_population}(:,1))
       xlabel('t');
   end
   title(['Population ',num2str(current_population),':  Total Population'],...
      'Color','blue','FontSize',10,'Units','normalized');
   set(axes_moment_handle,'Tag','Axes_population_moment');
   
   set(status_bar_handle,'String','Simulation Finished');
end
