%% LATEX TWICE TO GET THE REFERENCES RIGHT
%%\documentstyle[12pt]{article}
\documentstyle[12pt,epsf,psfig]{article}

\newcommand{\postscript}[2]
 {\setlength{\epsfxsize}{#2\hsize}
  \centerline{\epsfbox{#1}}}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newtheorem{theorem}{Theorem}
\newtheorem{itlemma}{Lemma}[section] %number by section (set in \em by default)
%%%same numbering counter as itlemmas:
\newtheorem{itproposition}[itlemma]{Proposition}
\newtheorem{itcorollary}[itlemma]{Corollary}
\newtheorem{itremark}[itlemma]{Remark}
\newtheorem{itnote}[itlemma]{Note}
\newtheorem{itobservation}[itlemma]{Observation}
\newtheorem{itdefinition}[itlemma]{Definition}
\newtheorem{itexercise}[itlemma]{Exercise}
\newtheorem{itexample}[itlemma]{Example}

\newenvironment{lemma}{\begin{itlemma}\rm}{\end{itlemma}} %no-italics
\newenvironment{remark}{\begin{itremark}\rm}{\end{itremark}} %no-italics
\newenvironment{observation}{\begin{itobservation}\rm}{\end{itobservation}} %no-italics
\newenvironment{corollary}{\begin{itcorollary}\rm}{\end{itcorollary}}
\newenvironment{proposition}{\begin{itproposition}\rm}{\end{itproposition}}
\newenvironment{definition}{\begin{itdefinition}\rm}{\end{itdefinition}}
\newenvironment{exercise}{\begin{itexercise}\rm}{\end{itexercise}}
\newenvironment{example}{\begin{itexample}\rm}{\end{itexample}}

\tolerance=1600 %allow some tolerance in extending after line.
\parskip=3pt plus 1pt minus 1pt
\overfullrule=0pt %no dark lines if overfull

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% This are for choices.  They take 4 or 6 args
%% (carriage returns allowed in between}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newcommand{\twoif}[4]{
\left\{ \begin{array}{ll}#1&#2\\#3&#4\end{array}\right.
}
\newcommand{\threeif}[6]{
\left\{ \begin{array}{ll}#1&#2\\#3&#4\\#5&#6\end{array}\right.
}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\def\partmlab#1#2#3#4#5#6#7#8{\bordermatrix{&#1&#2\cr#3&#5&#6\cr
                                        \noalign{\medskip}#4&#7&#8\cr}}
%% For 2by2 partitioned matrices; arguments are: two labels of top, side,
%% 4 entries (row by row).  If entries are themselves matrix or pmatrix,
%% better for readability of source to define first: \def\matrix1{...} etc.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\newcommand{\bigframebox}[1]{\framebox{\parbox{\smallerwidth}{#1}}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Warning: use \noindent before!
% AND IN TEXT, after begin doc, do:
% \newlength{\smallerwidth}  %for bigboxes
% \setlength{\smallerwidth}{\textwidth}
% \addtolength{\smallerwidth}{-8pt}        %MUST BE DONE AFTER begin doc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\def \beq {\begin{eqnarray}}
\def \eeq {\end{eqnarray}}
\def \beqn {\begin{eqnarray*}}
\def \eeqn {\end{eqnarray*}}

\newcommand{\rf}[1]{[\ref{#1}]}            %bibliographical references

\newcommand{\text}[1]{\hbox{\rm \ #1\ \/}}
\newcommand{\be}[1]{\begin{equation}\label{#1}}
\newcommand{\ee}{\end{equation}}
\newcommand{\bl}[1]{\begin{lemma}\label{#1}}
\newcommand{\ble}[1]{\begin{lemmaex}\label{#1}}
\newcommand{\br}[1]{\begin{remark}\label{#1}}
\newcommand{\bn}[1]{\begin{note}\label{#1}}
\newcommand{\bo}[1]{\begin{observation}\label{#1}}
\newcommand{\bt}[1]{\begin{theorem}\label{#1}}
\newcommand{\bd}[1]{\begin{definition}\label{#1}}
\newcommand{\bp}[1]{\begin{proposition}\label{#1}}
\newcommand{\bc}[1]{\begin{corollary}\label{#1}}
\newcommand{\bfact}[1]{\begin{fact}\label{#1}}
\newcommand{\ber}[1]{\begin{exercise}\label{#1}}
\newcommand{\bex}[1]{\begin{example}\label{#1}}
\newcommand{\bem}[1]{\begin{example}\label{#1}}    %Yes, 2 different ones...
\newcommand{\ec}{\mybox\end{corollary}}
\newcommand{\efact}{\mybox\end{fact}}
\newcommand{\eer}{\mybox\end{exercise}}
\newcommand{\eex}{\mybox\end{example}}
\newcommand{\eem}{\mybox\end{example}}
\newcommand{\el}{\mybox\end{lemma}}
\newcommand{\ele}{\mybox\end{lemmaex}}
\newcommand{\er}{\mybox\end{remark}}
\newcommand{\en}{\mybox\end{note}}
\newcommand{\eo}{\mybox\end{observation}}
\newcommand{\et}{\qed\end{theorem}}
\newcommand{\ed}{\mybox\end{definition}}
\newcommand{\ep}{\mybox\end{proposition}}
\newcommand{\epr}{\end{proof}}
\newcommand{\bpr}{\begin{proof}}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%  These are just for statements, no "qed" box is put in.
%%  (Also for when I have to finish something after a displayed equation
%%        or itemized.)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newcommand{\ecs}{\end{corollary}}
\newcommand{\eers}{\end{exercise}}
\newcommand{\eexs}{\end{example}}
\newcommand{\eems}{\end{example}}
\newcommand{\els}{\end{lemma}}
\newcommand{\eles}{\end{lemmaex}}
\newcommand{\ers}{\end{remark}}
\newcommand{\ens}{\end{note}}
\newcommand{\eos}{\end{observation}}
\newcommand{\ets}{\end{theorem}}
\newcommand{\eds}{\end{definition}}
\newcommand{\eps}{\end{proposition}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\newcommand{\qed}{\hfill \halmos} %put \qed at right margin
\newcommand{\mybox}{\hfill $\Box$} %put \qed at right margin (white square)

%%%\newcommand{\pmat}[2]{\left( \begin{array}{#1} #2 \end{array} \right) }
%%%use as e.g. \pmat{cc}{ a&b\\ a &c }  Note 'cc' for 2 Centered cols
%%%can also use \pmatrix{ a&b\cr a &c \cr} for plainTeX standard construct
%%%\setcounter{page}{27} should set the CURRENT page to 27

\def \R {\bb R} % reals
\def \C {\bb C} % complexes
\def \Z {\bb Z} % integers
\def \X {\bb X} % for this paper only!
\def \inter {{\rm int \,}}
\def \rank {{\rm rank \,}}
%%%\def \inf {{\rm inf \,}}  %better in Latex -- subscripts OK!
%%%\def \min {{\rm min \,}}
\def \im {{\rm im \,}}
\def \ker {{\rm ker \,}}
\def \norm #1{\Vert #1 \Vert}
\def \ip #1\by #2{\langle #1 , #2 \rangle}
\def \Sone {{\bf S}^1} % circle
\def \wigglearrow {\sim \mkern -6mu >}

\def \today {\ifcase \month \or January\or February\or  
  March\or April\or May\or June\or July\or August\or 
  September\or October\or November\or December\fi\   
  \number \day, \number \year}                       
\newcommand{\comment}[1]{}
\newcommand{\halmos}{\rule{1ex}{1.4ex}}
%%%alternative to \qed (see below)
\newenvironment{proof}{\noindent {\em Proof}.\ }{\hspace*{\fill}$\halmos$\medskip}

%%%\def\qed{\kern 1em \vrule height .5em width .5em depth 0em}
%%%(I now use new qed: in right margin)

\def\vbar{\mathchoice{\vrule height6.3ptdepth-.5ptwidth.8pt\kern-.8pt}
   {\vrule height6.3ptdepth-.5ptwidth.8pt\kern-.8pt}
   {\vrule height4.1ptdepth-.35ptwidth.6pt\kern-.6pt}
   {\vrule height3.1ptdepth-.25ptwidth.5pt\kern-.5pt}}
\def\fudge{\mathchoice{}{}{\mkern.5mu}{\mkern.8mu}}
\def\bbc#1#2{{\rm \mkern#2mu\vbar\mkern-#2mu#1}}
\def\bbb#1{{\rm I\mkern-3.5mu #1}}
\def\bba#1#2{{\rm #1\mkern-#2mu\fudge #1}}
\def\bb#1{{\count4=`#1 \advance\count4by-64 \ifcase\count4\or\bba A{11.5}\or
   \bbb B\or\bbc C{5}\or\bbb D\or\bbb E\or\bbb F \or\bbc G{5}\or\bbb H\or
   \bbb I\or\bbc J{3}\or\bbb K\or\bbb L \or\bbb M\or\bbb N\or\bbc O{5} \or
   \bbb P\or\bbc Q{5}\or\bbb R\or\bbc S{4.2}\or\bba T{10.5}\or\bbc U{5}\or
%%%   \bbb P\or\bbc Q{5}\or\bbb R\or\bba S{8}\or\bba T{10.5}\or\bbc U{5}\or
   \bba V{12}\or\bba W{16.5}\or\bba X{11}\or\bba Y{11.7}\or\bba Z{7.5}\fi}}
%
%
\renewcommand{\thefootnote}{\fnsymbol{footnote}}
%%% read as: theFOOTNOTE (i.e. the command  executed when doing a footnote) is
%%% now to be set as a FOOTNOTESYMBOL (type of number) value of counter
%%% FOOTNOTE (Instead of Roman).  See LaTeX manual pp.91-92 as well as 175.
\def\bi{\begin{itemize}}
\def\ei{\end{itemize}}
\def\ben{\begin{enumerate}}
\def\een{\end{enumerate}}
\def\st{\, | \,}

\newcommand{\myboxineq}{\eqno\Box} %put \qed at margin (white square) in eqn
\newcommand{\qedineq}{\eqno\halmos} %put \qed at margin in eqn

\topmargin -0.5in %default is 0.0, i.e. 1 inch from top.  Now it is 0.5
\textheight 9in
\oddsidemargin 0.25in %more reasonable margins!  (This is 1.25; can use <0 too)
\evensidemargin 0.25in %(same as above)
\textwidth 6.25in %(could use also -.5,-.5,7.5 for instance)
\tolerance=1600 %allow some tolerance in extending after line.
\parskip=3pt plus 1pt minus 1pt
\overfullrule=0pt %no dark lines if overfull

%%%% FROM HERE, macros special for paper
\def\mm{m}
\def\pp{p}
\def\proj{\pi _1} % projection on first group of coordinates
\def\x{u}  % vector argument
\def\SEC{(SEC)}
\def\SECo{(SEC$^0$)}
\def\INV{(INV)}
\def\Fp{{\cal F}_\pp}
\def\Fpm{{\cal F}_\pp^\mm}  %there is a 3,3 and a 2,1 in text

\def\Rm{\R^m}
\def\Rn{\R^n}
\newcommand{\ball}[2]{{\bb B}(#1,#2)}
\newcommand{\cball}[2]{\overline {\bb B} (#1,#2)}
\def\Fnm{{\cal F}_n^\mm}
\def\clos{{\rm clos\,}}
\def\inter{{\rm int\,}}

\def\CCF{compatible class of functions}
\begin{document}

\pagestyle{empty}
\begin{flushright}
{\sl SCCS report: 782}
\end{flushright}
\vfill
\begin{center}
\thicklines
\setlength{\unitlength}{0.6cm}
\begin{picture}(14.5,12)(-7.25,-6)

\multiput(-5.5,-1.45)(0,1.75){2}{\framebox(11,1.15){}}

\multiput(-1.45,-5.5)(1.75,0){2}{\framebox(1.15,11){}}

\multiput(-3.75,-2.35)(6.35,0){2}{\framebox(1.15,4.7){}}

\multiput(-2.35,-3.75)(0,6.35){2}{\framebox(4.7,1.15){}}

\multiput(-2.9,0)(2.9,-2.9){2}{\line(1,1){2.9}}

\multiput(-2.9,0)(2.9,2.9){2}{\line(1,-1){2.9}}

\end{picture}

\vspace{1cm}

\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\vspace{.05cm}
\hrule
\vspace{.25cm}
{\bf 
Mean Field Annealing and non-Euclidean Combinatorial Optimization
Problems: The case of Academic Course scheduling, (Part I: Algorithms,
Part II: Implementations)
}
\vspace{.25cm}
\hrule
\vspace{.05cm}
\hrule
\hrule
\hrule
\hrule
\hrule
\vspace{1cm}

Saleh Elmohamed\\
Geoffrey Fox\\
Paul Coddington\\
Dept. of Computer Science, and\\
Northeast Parallel Architectures Center,\\
Syracuse University\\
Syracuse, NY 13244\\
\{saleh, gcf, paulc\}@npac.syr.edu\\
\end{center}

\newpage
\pagestyle{empty}
\begin{flushright}
{\sl SCCS report: 782}
\end{flushright}
\vfill
\begin{center} 
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\vspace{.05cm}
\hrule
\vspace{.25cm}
{\bf 
Mean Field Annealing  and non-Euclidean Combinatorial Optimization
Problems: The case of Academic Course Scheduling, (Part I: Algorithms,
Part II: Implementations)  
}
\vspace{.25cm}
\hrule
\vspace{.05cm}
\hrule
\hrule
\hrule
\hrule
\vspace{1cm}

{\large
Saleh Elmohamed, Geoffrey Fox, and Paul Coddington\\}
\vspace{.5cm}
Northeast Parallel Architectures Center\\
Syracuse University, Syracuse, NY 13244\\
\vspace{.25cm}
\{saleh, gcf, paulc\}@npac.syr.edu\\
\end{center}

%\vspace{1cm}
%\centerline{\today}

\vspace{1.5cm}
\centerline{ABSTRACT}

In this study we have tackled the problem of academic course 
scheduling (or timetabling) at the university level using the method
of mean field annealing. We also extended the Potts normalization
conditions and equality constraints to accommodate our study case. 
Experimental analysis as well as summary of the results were included.

\vfill
\noindent{\bf Key words : feed-back neural nets, neurons, mean field annealing,  mean field theory, tanh functions, Potts systems, constraints, events, cost functions, phase transition.} 

% \vfill
% \centerline{\today}
% \vfill

% \newpage
% \pagestyle{empty}
% \tableofcontents

\newpage

\pagestyle{plain}
% \pagestyle{myheadings}
% \markright{\sl MFA Report}
\section{Introduction}
\setcounter{page}{1}

The prime objective of this study is to derive an approximate solution
to a non-geometric combinatorial optimization (CO)
problem, which is academic scheduling at the university level. 
It can be summarized as follows:
given data sets of classes and their days, enrollments, and
instructors; rooms and their capacities, types, and locations;
distances between buildings; priorities of each building for different
departments; and students and their class preferences; 
the problem is to construct a feasible class schedule satisfying all
the hard constraints and minimizing the  medium and soft constraints.
Hard constraints are space and time constraints that must be
satisfied, such as scheduling only one class at a time for any 
teacher, student, or classroom. 
Medium and soft constraints are student and teacher preferences that 
should be satisfied if possible.

The timetabling problem (TTP) is a high-dimensional, 
non-Euclidean, multi-constraint combinatorial optimization problem, 
and is consequently very difficult to solve.
This problem has been tackled by many researchers, 
mostly in the field of operations research.
A number of different heuristics have been tried on different 
instances of the problem, from high school to university course scheduling 
(see the reviews by de Werra~\cite{21} and Shaerf~\cite{50}) using 
method such as genetic algorithms \cite{37},  simulated annealing
\cite{2, 3, 39}, tabu search \cite{38},  neural networks \cite{32}, and
various combinations of two or more of these techniques. 
For small to medium size problems, some of these methods work well,
however no particular method has yet been shown to produce good results for 
large-scale, real-world problems. 
Very few of these studies have used real-world data, and 
none of these methods has been used on as large and complex a problem 
as we have addressed in this work. Also, we are not aware of any large
scale study that takes into account constraints due to student preferences.

We have used data for classes at Syracuse University. 
Currently this problem is handled by the university
scheduling department in a semi-automated fashion. A scheduling program 
is used to find a partial solution, and substantial manual effort is 
required to iterate towards a final solution. Also, when scheduling a
certain semester (e.g. fall 1996), a template of a previous semester
(e.g. fall 1995) is used as part of the input data.

We have applied a method consisting of  a feed-back artificial neural
system (ANS) composed of multi-state neurons using mean field theory
(MFT) equations to control the dynamics. This method is usually known
as mean field annealing (MFA). The MFT equations are used in place of the
computationally demanding stochastic updating procedures often used to
avoid getting stuck in local minima.

Much of the work in the neural networks field dealing with
optimization problems has centered on understanding and improving on
the Hopfield-Tank TSP algorithm \cite{32} (HT, for short). HT proposed
that CO problems could be solved using such analog networks. In
practice however, considerable difficulties have been encountered
carrying out this program for the TSP, leave alone problems such as
scheduling.

\section{The Timetabling Problem}

Timetabling is the assignment of time slots to a set of events,
subject to constraints on these assignments. 
The NP-complete {\it professors and classes} timetabling problem
\cite{2, 3, 27, 28} is a constraint satisfaction problem that can be
briefly stated as follows:

For a certain school with $N_p$ professors, $N_q$ classes, $N_x$
classrooms and lecture halls, and $N_s$ students, it is required to 
schedule $N_l$ professor-class pairs within a time
limit of $N_t$ time slots producing a legal schedule. A legal schedule
needs to be found such that no professor, class, or student is in more
than one place at a time, and no room is expected to accommodate more
than one lesson at a time or more students than its capacity. 
Figure \ref{relationship} sketches the approximate relationship
between the various entities of the timetabling problem.
\vspace{0.70cm}
\begin{figure}[htb]
\centerline{\psfig{figure=scheduling-net.eps,height=200pt,width=280pt}}
\vspace{0.20cm}
\caption{The relationship between various entities of TTP where P:
professors, C: classes, S: students, R: rooms, and I: time
periods. The direction of the arrow is the direction of the
assignment.}
\label{relationship}
\end{figure}

The constraints for this problem can be hard, medium or soft. 
The medium and soft constraints have an associated cost (or penalty), 
and if they are not satisfied, the goal is to minimize this cost. 
Soft constraints have a lower priority (and thus lower cost) than medium 
constraints.
The hard constraints must be satisfied, so their associated cost must be 
reduced to zero. A feasible schedule is one that satisfies all the hard
constraints.

\noindent
{\bf Hard constraints} are usually constraints that physically cannot be
violated. This includes events that must not overlap in time, such as:
\begin{itemize}
\item classes taught by the same professor,
\item classes held in the same room,
\item a class and a recitation or a lab of the same class.
\end{itemize}
%
Another examples are space or room constraints: 
\begin{itemize}
\item A class cannot be assigned to a
particular room unless the capacity of the room  
is greater than or equal to the class enrollment.
\item
Some classes, such as laboratories, require a certain type of room.
\end{itemize}

\noindent
{\bf Medium constraints} are usually considered to those constraints that fall 
into the gray area between the hard and soft constraints \cite{49}. 
In our implementation, we define medium constraints to be constraints such as 
time and space conflicts which, like hard constraints, cannot physically be
violated (for example, it is not possible for one person to be in two
different classes at the same time). However we consider these constraints 
to be medium rather than hard if they can be avoided by making adjustments
to the specification of the problem. The primary example is student 
preferences. We cannot expect to be able to satisfy all student class 
preferences, in some cases, certain students will have to adjust their 
preferences since certain classes will clash, or will be oversubscribed.

Medium constraints have a high penalty attached to them, although not as
high as that associated with the hard constraints. In the final
schedule the penalty of these constraints should be minimized and 
preferably reduced to zero. Some examples of medium constraints are:
\begin{itemize}
\item Avoid time conflicts for classes with students in common.
\item Eligibility criteria for the class must be met.  
\item Do not enroll athletes in classes that conflict with their sport
      practice time (of course, depending on the sport).
\end{itemize}


\noindent
{\bf Soft constraints} are preferences that do not deal with time conflicts, 
and have a lower penalty (or cost) associated with them. We aim to 
minimize the cost, but do not expect to be able to reduce it to zero.
Some examples are:
\begin{itemize}
\item For each student, balance the three-day ({\it Mon, Wed, Fri}) 
      as well as the two-day ({\it Tue, Thu}) schedules.
\item Balance or spread out the lectures over the week.
\item Classes may request contiguous time slots.
\item Balance enrollment in multi-section classes.
\item Lunch and other break times may be specified.
\item Professors may request periods in which their classes are not taught.
\item Professors may have preferences for specific rooms or types of rooms.
\item Minimize the distance between the room where the class is assigned and 
      the building housing its home department. 
\end{itemize}

Some soft constraints may have higher priority (and thus higher cost) than
others. For example, preferences involving teachers will have higher priority 
than the preferences of students.

The {\bf cost function} measures the quality of the current schedule and
generally involves the weighted sum of penalties associated with
different types of constraint violations. The aim of the optimization 
technique is to minimize the cost function.

% \section{The Constraint Satisfaction Problem}
%
% We now briefly describe the general constraint satisfaction problem
% which has been widely studied in the Artificial Intelligence community
% \cite{9, 10, 11}.
% 
% {\it An instance} of a constraint satisfaction problem {\it consists
% of}:
% 
% \begin{itemize}
% \item{\it A finite set of variables, $N_i$, where $i = 1, 2,\ldots,n\;$}.
% \item{\it A domain of D of possible values of $N_i$}$\;$.
% \item{\it A set of constraints} $\{C(S_1), C(S_2),\ldots, C(S_c)\}\;$.
%
% {\it Each $S_i$ is an ordered subset of the variables, and each
% constraint $C(S_i)$ is a set of tuples indicating the mutually
% consistent values for the variables in $S_i$}$\;$.
% \end{itemize}
%
% The length of the tuples in a given constraint will be called the
% {\it arity} of that constraint\footnote{Strictly speaking we also need
% to associate an arity with an empty set of tuples, but we will neglect
% this special case, in order to simplify the presentation.}. In
% particular, unary constraints specify the allowed values for a single
% variable, and binary constraints specify the allowed combinations of
% values for a pair of variables.
%
% A {\it solution} to a constraint satisfaction problem is an assignment
% of values to the variables which is consistent with all of the
% constraints. Deciding whether or not a given problem instance has a
% solution is NP-complete in general \cite{11}, even when the
% constraints are restricted to binary constraints.

\section{Brief Notes on the Complexity of Timetabling}
Timetabling belongs to the class of constraint satisfaction problems
and it is essentially an example of a resource-constrained scheduling
problem, which is an NP-hard \cite{12, 58}. Here, resources
are physical entities, for example: students, rooms, overhead
projectors, and so on. Also time periods can be thought of as
resources.

In general, timetabling problems are always NP-complete
\cite{25}, which means there exists no known polynomially bound
algorithm for solving them optimally. As a result, these problems are
often solved by means of {\it heuristics} -- solution procedures that
focus on finding a feasible schedule of ``good'' (as opposed to
optimal) quality within an acceptable amount of time.

Traditionally, timetabling has been approached by means of linear
programming (LP) with binary variables plus some heuristics. For
example, if $i$ identifies a teacher, $j$ identifies a time interval
and $k$ identifies a class then the binary variable $X_{ijk}$ $=$ $1$
if teacher $i$ has class $k$ at time-interval $j$; $X_{ijk}$ $=$ $0$ ,
otherwise. Suppose we have a set of data as small as 20
teachers, 30 time intervals, and 10 classes. If we use the LP approach
to obtain a feasible schedule out of this data, then we need to deal
with over 6000 variables to represent the problem, and 
$2^{6000}$ possible states. The problem quickly becomes intractable as
the number of variables increases.

For the much simpler problem of {\it class allocation}, for which there
is no timetabling involved (i.e. no allocation of classes to time 
periods), the computational load ($n_c$) from tackling this problem
using the LP method scales as follows \cite{5, 28}: 

\be{complexity}
n_c \propto N_c N_x N_t\;\;,
\ee

where $N_c$ is the number of room categories available
(i.e. classroom, auditorium, etc.).

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Mean-Field Annealing (MFA)}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
One of the potential drawbacks of using methods such as simulated
annealing for hard  optimization problems is that finding a good
solution can often take an  unacceptably long time. Mean-field
annealing attempts to avoid this problem by using a
deterministic approximation to simulated annealing, by attempting to
average over the statistics of the annealing process.  The result is
improved execution speed at the expense of solution quality.  Although
not strictly a continuous descent technique, MFA is closely related to
the Hopfield neural network~\cite{32}. 

MFA has been successfully applied to high-school class  scheduling
\cite{28}. For scheduling, it is advantageous to use a  Potts neural
encoding to specify discrete neural variables (or neurons)  for the
problem. This is defined in its simplest form as a mapping  of events
onto space-time slots, for example events defined by a professor-class
pair ($p,q$) are mapped onto  space-time slots ($x,t$). Potts neurons
$s_{pq,xt}$ are defined to be 1 if event ($p,q$) takes place in the
space-time slot ($x,t$), and 0 otherwise. In this way, the constraints
involved can be embedded in the neural net in  terms of a Potts
normalization condition such as $\sum_{x,t} S_{pq,xt} = 1$.

For a full derivation of the mean-field annealing algorithm from its
roots in statistical physics, see Peterson {\it et al.} \cite{26}.
Here, we will summarize the algorithm as well as Hopfield's model.

\subsection{Some Background on the Hopfield Model}

The network is built by connecting a large number of neurons to
each other to form a dense, parallel array and, in general, any
$i^{th}$ neuron in the net is described by two variables: its current
state (input) $u_i$,  and its output $v_i$.
The output is related to the input by a monotonically increasing
transfer function: $v_i = g_{\lambda}(u_i)$. 
For large $\lambda$, $g_{\lambda}(u_i) = {1 \over {1 + e^{-2\lambda u_i}}}$
approximates a step function, forcing $v_i$ to
$[0,1]\;\Longrightarrow\;$ (the output of the net lies on the unit
hypercube). 

\subsubsection{Operation of the Network}

\vspace{0.1cm}
\begin{figure}[htb]
\begin{center}
\begin{picture}(320,330)
\psfig{figure=hopfield_net.eps,height=320pt,width=330pt}
\end{picture}
\end{center}
\caption{A schematic view of the Hopfield network.}
\label{Hop-model}
\end{figure}

The network (see Figure \ref{Hop-model}) operates as follows:
Its dynamics are governed by a set of nonlinear coupled differential
equations of the form: 

\be{FODE.2}
{{du_i} \over {dt}} = - u_i + \sum_{j=1}^{n} W_{ij} v_j + I_{i}^{b}\;.
\ee

or the general form of it:

\be{FODE.3}
{{du} \over {dt}} = - \eta u + W v + I^b\;.
\ee

where $W_{ij}$ is the connection strength between neuron $i$ and neuron
$j$, $I_i^b$ is an input bias to neuron $i$, and the $\eta$ term
introduces an element of decay into the network's dynamics, which can
be used to ensure that ${\bf u}$ variables remain bounded in magnitude.

\subsubsection{Updating of the States of the Neurons}

To make sure that stability of the of the net is guaranteed, $u_i$ is
updated continuously by evaluating ${d \over dt}u_i$ as a function of
the total input to the $i^{th}$ neuron. 

\be{State.1}
 \dot{u} = 
-\eta u +  W v +  I^{b}\;,
\ee

where 

\be{State.2}
v = g(u)\;.
\ee
Equations \ref{State.1} and \ref{State.2} govern the dynamics of the
network.

\subsubsection{Stability of the Net}

Given that $W$ is symmetric, Hopfield proved stability for the
network by showing that with the above dynamic equations
(\ref{State.1} and \ref{State.2}) the network output $v$ evolves
so as to minimize the {\it Liapunov function}:

\be{Liapunov}
E^{liap}(v) = 
\left( -{1 \over 2} \sum_{ij} W_{ij} v_i v_j -\sum_i I_{i}^b v_i \right) 
+ \eta{\displaystyle \sum_i}\, {\displaystyle \int_{0}^{v_i}}\, g^{-1}(V)\,dV\;.
\ee

$v = (v_1 \ldots v_n)$ is denoting the net state, therefore we
have an energy function defined throughout $v \in [0,1]^n$ and
dynamics which guarantee convergence to minima. For a suitable $W$,
the minima of the first term lie at the hypercube vertices. The second
term is minimized at the hypercube center but has negligible impact at
low $\eta$. 

\subsubsection{Updating the States of the Net}

Updating of the states (or the output of the neurons) is carried out
either asynchronously or synchronously.  From a simulation point of
view, we found that synchronous updating a bit more efficient. From a
theoretical point of view, it makes not much of a difference on the
operation of the network as to which updating mode is used.

\subsubsection{Mapping the Problem onto the Network}

First, a representation of the problem is needed in which the feasible
solutions lie at hypercube vertices.  Each hypercube vertex is
regarded as a configuration of the problem, with an associated energy
(or cost) function. The partition function is then formed by summing
the Boltzmann factors of the admissible configurations. The main part
to derive the solutions is a quadratic energy function (to be mapped
onto the first term of \ref{Liapunov}) for which the minima correspond
to solutions and the depth of each minimum reflects the solution quality. 

\subsubsection{Shortcomings of the Hopfield Model}

The network's partition function sums over a vast number of
configurations (e.g. for the TSP, all of the ${2^{N}}^2$
configurations are admissible). These configurations, again, for the
TSP, are nothing like valid tours, and even though these offenders
have small Boltzmann weights their large number inevitably affects the
thermal average quantities.

\section{The Potts Approach to MFA for Timetabling}

Clearly the original HT algorithm had difficulty constraining the
network into valid solutions. One way forward is therefore to lessen
the computational burden being placed on the network. This can be
done by constraining each city (in the TSP case) to be ``on'' only once,
by enforcing $\sum_a V_{ia} = 1\; \forall i$, ($V_{ia} = 1/0$ means
that the city $i$ is (is not) the $a^{th}$ city visited in the
tour), rather than relying on an energy penalty term to try and do
this (as in HT). This idea was first utilized and analyzed for the TSP
case by, among others, Peterson {\it et al.} \cite{26}. It greatly
enhances the degree of convergence to valid solutions.
Later, Peterson and colleagues somewhat extended the technique and
applied it to a case of high school scheduling \cite{27, 28}.

In \cite{26} the normalization occurs in
the context of a {\it mean field annealing} approach. This generally
yields the same final network equations as the ``neuronal circuit''
approach of \cite{33}, but the exposition is a little bit clearer, as
it is laid out  in a statistical mechanics framework\footnote{This
equivalence was also noted in \cite{32}. P. Simic and Martin Simmen
brought this point to our attention.}. The idea is to regard each
hypercube vertex (see \cite{33}) as a configuration, with an
associated energy (or cost) ${\it E}$. The partition function is then 
formed by summing the Boltzmann factors of the  admissible
configurations. For example, in the HT algorithm, for the TSP, as
previously noted all of the $2^{N^2}$ configurations are admissible,
whereas if each city is restricted to being visited only once, then
only $N^N$ vertices are admissible. After taking a mean field
approximation, saddle point equations are derived (see next sections),
the solutions of which to pick out the dominant states of the network
at the current temperature $T$.

This statistical mechanics characterization gives a somewhat clearer
picture of why the normalized model (which corresponds
to a type of Potts model \cite{36} in physics, in that each
city represented by a spin which can be in only one of $N$ states)
ought to perform better than HT model. As Peterson {\it et al.} \cite{26}
indicated, the HT algorithm's partition function sums over a vast
number of configurations which are nothing like valid tours in the TSP,
and even  though these offenders have small Boltzmann weights their
large number inevitably affects the thermal average quantities.
The other advantage that draws us to this derivation is that, by
introducing a temperature, the concept of annealing, i.e. reducing
$T$ during a run, can be justifiably employed.

The above two paragraphs outline more or less the reasons
behind our choice to use MFA.

\subsection{Constraint Formulations in the Potts Model}

In this section we define some key terms and concepts used 
in this paper and following the notation of Peterson and colleagues on
this subject.

Let ($N_p$) to be a number of professors that give ($N_q$) classes in 
($N_x$) class-rooms at ($N_t$) time slots. Event-wise, $(p,q)$ event
-- denoting a professor $p$ giving a class $q$ -- takes place in the
space-time slot $(x,t)$ -- denoting classroom $x$ at time $t$.
Let $S_{pq;xt}$ denote the Potts (or multi-state) neuron, where
\be{form-1}
S_{pq;xt} =
\twoif {1} {\text{if event $(p,q)$ takes place in $(x,t)$ ,}}
          {0} {\text{otherwise.}}
\ee

Now, we define a number of {\bf hard  constraints} as follows: It is
known that an event $(p,q)$ should occupy precisely one space-time
slot $(x,t)$ so its associated constraint can be embedded in the
neural net in terms of the Potts normalization condition: 
\be{form-2}
\sum_{x,t} S_{pq,xt} = 1
\ee
Different events $(p_1,q_1)$ and $(p_2,q_2)$ should not occupy the same
space-time slot $(x,t)$.
\be{form-3}
E_{XT} =
{1 \over 2}\sum_{x,t}\sum_{p_1,q_1}\sum_{p_2,q_2}S_{p_1q_1;xt}S_{p_2q_2;xt} = 
{1 \over 2}\sum_{x,t}{\left(\sum_{pq}S_{pq;xt}\right)}^2
\ee
A teacher $p$ should have at most one class at a time.
\be{form-4}
E_{PT} =
{1 \over 2}\sum_{p,t}\sum_{q_1,x_1}\sum_{q_2,x_2}S_{pq_1;x_1t}S_{pq_2;x_2t} = 
{1 \over 2}\sum_{p,t}{\left(\sum_{qx}S_{pq;xt}\right)}^2
\ee
A class $q$ should have at most one teacher at a time.
\be{form-5}
E_{QT} =
{1 \over 2}\sum_{q,t}\sum_{p_1,x_1}\sum_{p_2,x_2}S_{p_1q;x_1t}S_{p_2q;x_2t} = 
{1 \over 2}\sum_{q,t}{\left(\sum_{px}S_{pq;xt}\right)}^2
\ee

In \cite{32, 33},  Hopfield introduced a model with neurons having a graded
response. So given $I_i$ a fixed input bias to neuron $x_i$, $W_{ij}$ the weight
of the connection between neurons $(x_i,x_j)$, and parameter 
$T$ to denote the temperature, the behavior of $x_i$ can be stated as:
 \be{form-JH}
 x_i =
 {1 \over 2} + {1 \over 2} \cdot 
 tanh\left(({\sum_{i \neq j}W_{ij}x_j - I_i}) \times {1 \over T} \right)
 \ee
When $T \longrightarrow 0$, function ({\ref{form-JH})'s
response takes on a step-wise form. These models of \cite{33}
are known to be a bit sensitive to parameter tuning for real world
applications, especially those of the type that involves ``finding''
or ``searching for'' $1$ out-of $N$ type problems.  

At a finite ${T}$, $x_i$ (of eq. \ref{form-JH})
can be replaced by an MFT variable $\upsilon_{ia}$, then
differentiating to obtain the Potts MFT equations:

\be{form-5.1}
\upsilon_{ia} =
- {1 \over {T}} \cdot {{\partial E}[{\vec\nu}] \over {\partial \nu_{ia}}}
\ee
\be{form-5.2}
\nu_{ia} =
{e^{\upsilon_{ia}} \over {\sum_b e^{\upsilon_{ib}}}}
\ee

Which enforces

\be{form-5.3}
\sum_{a} \nu_{ia} = 1\,\,\,\, \forall i\, .
\ee
Where $a$ denotes the vector components or all the elements of the group.
Also, we could state that $\nu_{ia}$ $\approx$ Probability(neuron(i) $\in$ state(a)).
Equations (\ref{form-5.1}, \ref{form-5.2}) are computed in an iterative manner.

In \cite{26, 18}, it was stated that the performance of MFA on certain
problems can be improved using  {\it neuron normalization}. This idea
was originally inspired by an analogy with Potts glasses in
statistical physics. {\it Neuron normalization} modifies the MFA
algorithm so that constraints of the form  $\sum_{a \in \text{States}}
\nu_{ia} = 1$ are explicitly enforced, without the need for a penalty
term. It was found by Peterson and other to be particularly useful for
mapping the TSP and graph partitioning (GP) type problems.

Now, to carry out the minimization of the energy functions stated
above, we need to convert the neural variables
$\{0,1\}$ of $S_i$ to their corresponding {\it mean field variables}
at $T$. These correspondences are derived (see \cite{28}) to be: 

\be{form-6.1}
V_{pq;xt} \approx <S_{pq;xt}>_{T}\,\,;
\ee
which is the thermal average of the original binary spines, with the following 
multi-state $MFT$ equations (see (\ref{form-5.1}, \ref{form-5.2}) for 
a general format):

\be{form-6}
U_{pq;xt} =
- {1 \over {T}} \cdot {{\partial E} \over {\partial V_{pq;xt}}}
\ee
\be{form-7}
V_{pq;xt} =
{e^{U_{pq;xt}} \over {\sum_{x^{'}t^{'}}e^{U_{pq;x^{'}t^{'}}}}}
\ee

The $MFT$ equations of $\upsilon$ and $\nu$ are
solved {\it iteratively} either synchronously or asynchronously, under
${T}$. Also, these equations give rise to 2 phase  transitions; one in
$x$ (space) and the other in $t$ (time). This leads to factorization
of the system by replacing $S_{pq;xt}$ by {\bf $X$-neurons}
$S^{(X)}_{pq;x}$ and {\bf $T$-neurons} $S^{(T)}_{pq;t}$ :

\be{form-8}
S_{pq;xt} =
S^{(X)}_{pq;x} \cdot S^{(T)}_{pq;t}
\ee
With separate Potts conditions:
\be{form-9}
\sum_{x}S^{(X)}_{pq;x} = 1
\ee
\be{form-10}
\sum_{t}S^{(T)}_{pq;t} = 1
\ee
For {\bf $X$-neurons} and {\bf $T$-neurons}, respectively. 
Rather than dealing with $p$ and $q$ directly, an independent variable
($i$) -acting as an event index- to which $p$ and $q$ are attributes
$p(i)$ and $q(i)$. The value of these exits in a look-up table which
contains all the relevant information to process each event ($i$).
Now, the above (\ref{form-9}, \ref{form-10}) conditions can be 
re-stated as follows:

\be{form-9.1}
\sum_{x}S^{(X)}_{i;x} = 1
\ee
\be{form-10.1}
\sum_{t}S^{(T)}_{i;t} = 1
\ee

During annealing, as we go from a higher temperature to a lower one, 
a {\it phase transition} is passed at (${T} = {T}_c$), where 
${T}_c$\footnote{See pp 29-31, \cite{1} for an approximate
value in the case of Ising ferromagnets.} is the {\it critical temperature}. 
The value of ${T}_c$ set the scale or the range of the value of
${T}$, and it is estimated by expanding \ref{form-7} around some
trivial fixed-point such as

\be{trivial-form}
V^{(0)}_{pq;xt} = {1 \over {N_x N_t}} \;\text{, at ${T}$ = 0} .
\ee

Also ${T}_c$ depends on the self-coupling coefficient $\beta$, and on
the method of the update.  This parameter $\beta$ acts as a Lagrangian
multiplier for the given energy constraint(s). See section on the
generic Potts algorithm and \cite{26}.

During the revision process of an existing schedule, sometimes it is
desirable to have a certain cluster of events to stick  together in
time (or within some specified period of time). This is generally
referred to as {\it relative clamping} and amounts to making the
following notational replacement: 
\be{form-12}
i \longrightarrow (j,k) \;\; ;
\ee
where $(j)$ denotes the cluster and $(k)$ an event within $(j)$. In this
situation, one has a common $t$-neuron $(S^{(T)}_{j;t})$ for cluster
$j$ but distinct $x$-neurons $(S^{(X)}_{jk;x})$ for the individual events.

To deal with events that take place over two consecutive hours, what is
called effective $time$-neurons $\tilde{S}^{(T)}_{i;t}$ are introduced:
\be{form-11}
\tilde{S}^{(T)}_{i;t} =
\sum_{k=0}^{g_i - 1} S^{(T)}_{i;t-k} \;\;{\text{;}}
\ee
where $g_i = 1$ for single hours and $= 2$ for double hours. For those 
events that occur between one and two hours, $g_i$ can be adjusted 
accordingly to handle these cases.
Using these effective $time$-neurons, constraints (\ref{form-3},
\ref{form-4}, \ref{form-5}) can be put into the following forms: 

\be{form-11.1}
 E_{XT} =
 {1 \over 2}\sum_{x,t}\sum_{i \neq i^{'}}
 S^{(X)}_{i;x}S^{(X)}_{i^{'};x} \tilde{S}^{(T)}_{i;t}
 \tilde{S}^{(T)}_{i^{'};t}
\ee

\be{form-11.2}
 E_{PT} =
 {1 \over 2}\sum_{p,t} \sum_{i \neq i^{'}} \delta_{p(i),p(i^{'})}
 \tilde{S}^{(T)}_{i;t} \tilde{S}^{(T)}_{i^{'};t}
\ee

\be{form-11.3}
 E_{QT} =
 {1 \over 2}\sum_{q,t} \sum_{i \neq i^{'}} \delta_{q(i),q(i^{'})}
 \tilde{S}^{(T)}_{i;t} \tilde{S}^{(T)}_{i^{'};t}
\ee

The Kronecker $\delta_{q(i),q(i')}$ ensures that only
events $i$ and $i'$ with identical ``primordial'' classes (or point of
origination) are,  for example, summed over. Here, $q(i)$ and $q(i')$
refer to classes of the same type or 2 sections derived from the same
original class.
\be{Kronecker}
\delta_{r,s} =
\twoif {1} {\text{if} r = s}
       {0} {\text{otherwise.}}
\ee


One important soft constraint that need to be handled is {\it group formation}.
Briefly, it can be stated as follows: for some optional subjects the
classes are broken up into option groups temporarily forming new 
classes. In order to account for this, we need a formalism that allows
for the breaking up of these primordial classes and subsequent 
recombinations into quasi-classes. To carry out such an extension one
might have contributing pairs of events where the primordial 
classes are different. Hence one should replace $\delta$ with
a more structured overlap matrix $\Gamma$ which is given by a look-up
table.
\be{form-11.4}
\delta_{q(i),q(i^{'})} \longrightarrow \Gamma_{q(i),q(i^{'})} \; .
\ee
Using relative clamping and form \ref{form-11.4}, we need to make 
sure that ($\forall k$) events within a given ($j$) cluster $p(j,k)
\neq 0$. This says that there can be no ``uncovered'' class period or
a class event (that is, periods when no professor has been assigned to
an already scheduled class). For this purpose we introduce
$\tilde{\delta}_{j,k}$ according to:
\be{form-3.9.3}
\tilde{\delta}_{j,k} = \delta_{j,k}(1 - \delta_{j,0}\delta_{k,0}) \; .
\ee
If $j = 0$ then $\forall k$ , $k = 0$ .

Using the abovementioned concepts, the collision terms of equations 
(\ref{form-11.1}, \ref{form-11.2}, and \ref{form-11.3})
are extended to ensure that all the hard constraints are satisfied (or
the solutions obtained are legal). The extension takes the following
forms:

\be{form-3.9.4.1.1}
 E_{XT}{\text{(form 1)}} =
 {1 \over 2}\sum_{jk \neq j^{'}k^{'}}\sum_{x}
 S^{(X)}_{jk;x}S^{(X)}_{j^{'}k^{'};x}
\sum_{t}\tilde{S}^{(T)}_{j;t}\tilde{S}^{(T)}_{j^{'};t}
\ee

\be{form-3.9.4.1.2}
 E_{XT}{\text{(form 2)}} =
 {1 \over 2}\sum_{j \neq j^{'}}\sum_{k k^{'}}\sum_{x}
 S^{(X)}_{jk;x}S^{(X)}_{j^{'}k^{'};x}
 \sum_{t}\tilde{S}^{(T)}_{j;t}\tilde{S}^{(T)}_{j^{'};t}
 + {1 \over 2}\sum_{j}g_j \sum_{k \neq k^{'}} \sum_{x} 
 S^{(X)}_{jk;x}S^{(X)}_{jk^{'};x}
\ee

\be{form-3.9.4.2.1}
 E_{PT}{\text{(form 1)}} =
 {1 \over 2}
\sum_{p,t} \left[ \sum_{jk} \tilde{\delta}_{p,p(j,k)}
 \tilde{S}^{(T)}_{j;t} \right]
\left[ \sum_{j^{'}\neq j}\sum_{k^{'}}
 \tilde{\delta}_{p, p(j^{'},k^{'})} 
 \tilde{S}^{(T)}_{j^{'};t} \right]
\ee

\be{form-3.9.4.2.2}
 E_{PT}{\text{(form 2)}} =
        {1 \over 2}\sum_{j \neq j^{'}}\sum_{k k^{'}} 
        \tilde{\delta}_{p(j,k), p(j^{'},k^{'})}
        \sum_{t}\tilde{S}^{(T)}_{j;t} \tilde{S}^{(T)}_{j^{'};t}
\ee

\be{form-3.9.4.3.1}
 E_{QT}{\text{(form 1)}} =
 {1 \over 2} 
\sum_{qt} \left[ \sum_{jk} \delta_{q,q(j,k)} \tilde{S}^{(T)}_{j;t} \right]
\left[ \sum_{j^{'} \neq j} \sum_{k^{'}} \Gamma_{q,q(j^{'}, k^{'})}
        \tilde{S}^{(T)}_{j^{'};t} \right]
\ee

\be{form-3.9.4.3.2}
E_{QT}{\text{(form 2)}} =
        {1 \over 2}
        \sum_{j \neq j^{'}}\sum_{kk^{'}} 
        \Gamma_{q(j,k), q(j^{'},k{'})} \sum_{t} \tilde{S}^{(T)}_{j;t}
        \tilde{S}^{(T)}_{j^{'};t}
\ee

Note: Each (form 2) above is a derivation from its corresponding (form 1).

The factorization idea shown in equations 
(\ref{form-3.9.4.1.1}, \ref{form-3.9.4.1.2}, \ref{form-3.9.4.2.1},
\ref{form-3.9.4.2.2}, \ref{form-3.9.4.3.1}, and \ref{form-3.9.4.3.2}) 
can be extended to {\it mean field} variables $V^{(X)}_{jk;x}$ and
$V^{(T)}_{j;t}$, corresponding to neurons $(S^{(X)}_{jk;x})$ and
$(S^{(T)}_{j;t})$, respectively. $V^{(X)}_{jk;x}$ and $V^{(T)}_{j;t}$
can be interpreted as being the probabilities that the corresponding
events occur at given $x$-value and $t$-value with normalization given
by the equations (\ref{form-9.1}, \ref{form-10.1}) .

Substituting $V$ for $S$ in equations (\ref{form-9.1}, \ref{form-10.1}), 
the {\it mean field} equations at a finite  ${T}$, are given in terms
of the {\it local fields}  $U^{(X)}_{jk;x}$ and $U^{(T)}_{j;t}$. Now,
the {\it mean field} dynamics proceed over the following local and
mean field variables:

\be{form-4.1.3}
U^{(X)}_{jk;x} =
- {1 \over {T}} \cdot {{\partial E} \over {\partial V^{(X)}_{jk;x}}}
\ee

\be{form-4.1.4}
U^{(T)}_{j;x} =
- {1 \over {T}} \cdot {{\partial E} \over {\partial V^{(T)}_{j;t}}}
\ee

\be{form-4.1.5}
V^{(X)}_{jk;x} =
{e^{U^{(X)}_{jk;x}} \over {\sum_{x^{'}} e^{U^{(X)}_{jk;x^{'}}}}}
\ee

\be{form-4.1.6}
V^{(T)}_{j;t} =
{e^{U^{(T)}_{j;t}} \over {\sum_{t^{'}} e^{U^{(T)}_{j;t^{'}}}}}
\ee

The iterations of the {\it MF} dynamics over equations 
({\ref{form-4.1.3}} to {\ref{form-4.1.6}})
can either be {\bf synchronous} or {\bf serial} updating.
Synchronous updating involves a single sweep consists of first
computing all the local fields, and then updating all the
neurons. That means first perform equations
({\ref{form-4.1.3}},{\ref{form-4.1.4}}) then perform equations
({\ref{form-4.1.5}},{\ref{form-4.1.6}}). While serial updating is done
as follows: for each neuron $S^{(X)}_{jk;x}$ and/or $S^{(T)}_{j;t}$,
the {\it local field} is computed immediately before the corresponding
{\it neuron state} is updated. That means either do equation
({\ref{form-4.1.3}}) followed by ({\ref{form-4.1.5}}) and/or do
equation ({\ref{form-4.1.4}}) followed by ({\ref{form-4.1.6}}). 

Next, we need to formalize two important {\bf soft constraints} encountered
when scheduling classes.

{\em Spreading}: different classes in a particular area (i.e. computer science) 
should be spread over the week in a (MWF) and (TuTh) type of spread.
First, assume we have access to a subject attribute $s(j,k)$. Second, we
define an effective ``day-neuron'' as follows 
($\tilde{S}^{(T)}_{j;dh}$ is reduced to a $D$-type neuron by summing over 
all hours of the day):
\be{day-neuron}
S^{(D)}_{i;d} = \sum_h \tilde{S}^{(T)}_{j;dh} = g_j \cdot S^{(T)}_{j;dh}
\ee

So the penalty term $E_{QSD}$ that spreads classes $Q$ in a particular subject
$S$ over different week-days $D$ :

\be{spread-1}
E_{QSD} = 
{1 \over 2} \sum_{qsd}{\left[\sum_{jk}\sum_{h} \delta_{q,q(j,k)} \delta_{s,s(j,k)}
\tilde{S}^{(T)}_{jk;dh}\right]}^2
\ee

Refining the above a bit more by considering 2 different clusters $(j)$ and $(\hat{j})$ 
of classes over
$D$ : ($q$ and $s$ are not used in the summation, since no change occurs to them)

\be{spread-2}
E_{QSD} =
{1 \over 2} \sum_{d}\sum_{j \neq \hat{j}}\sum_{kh}\sum_{\hat{k}\hat{h}}
\delta_{q(j,k),q(\hat{j},\hat{k})}\delta_{s(j,k), s(\hat{j},\hat{k})}
\hat{S}^{(T)}_{j;dh}\hat{S}^{(T)}_{j;d\hat{h}}
\ee

Now for some day $(d)$ of the week:

\be{spread-3}
E_{QSD} =
{1 \over 2} \sum_{j \neq \hat{j}} \sum_{k\hat{k}} 
\delta_{q(j,k),q(\hat{j},\hat{k})}\delta_{s(j,k), s(\hat{j},\hat{k})}
\sum_d S^{(D)}_{j,d} S^{(D)}_{\hat{j},d}
\ee

{\em Gluing}: The schedule should have as few ``holes'' as possible; that is
lessons should be glued together. We would like to point out that this
constraint operates directly opposite to {\em spreading}, as a result
when {\em spreading} is present we have not dealt with {\em gluing}
and vice versa.

The penalty term that rewards situations where $(d,h)$-events are glued to
$(d,h-1)$-events can be defined as follows:

\be{glueing-1}
E_{QDH} = -\alpha \sum_{q}\sum_{d,h} \left[\sum_{j,k}\delta_{q,q(j,k)}\tilde{S}^{(T)}_{j;dh}\right]
\left[\sum_{\hat{j} \neq j}\Gamma_{q,q(\hat{j},\hat{k})} \tilde{S}^{(T)}_{j;d(h-1)}\right] 
\ee

The energy term of (\ref{glueing-1}) can be summarized to say that:
{\it if cluster ($j$) with events each symbolized by ($k$) occurs on
  hour ($h$) of day ($d$) then it would be  quite desirable to have
  another event ($\hat{k}$) from a different cluster ($\hat{j}$) to
  occur on hour ($h-1$) of that day.}

Next, since ($q$) is fixed we could re-write the above equation as follows:

\be{glueing-2}
E_{QDH} = -\alpha\sum_{\hat{j} \neq j}\sum_{k,\hat{k}}\Gamma_{q(j,k),
q(\hat{j},\hat{k})} \sum_{d,h>1}\tilde{S}^{(T)}_{j;dh}\tilde{S}^{(T)}_{\hat{j},d(h-1)}
\ee
The presence of the ($-$) sign in front of (\ref{glueing-1}
and \ref{glueing-2}) indicates that such energy functions are to be minimized. 
Also, the ($\alpha$) in front of the the above forms governs the strength of 
this reward relative to the energies of equations
(\ref{form-3},\ref{form-4},\ref{form-5}). This ($\alpha$) takes on a
value between zero and $1$.

Finding {\em legal} solutions with no collisions corresponds to minimizing
the {\em hard} energy:
\be{Hard}
E_{hard} = combine (E_{XT}, E_{PT}, E_{QT})
\ee
to zero. The {\em soft energy} is given by:
\be{Soft}
E_{soft} = combine (E_{QSD}, E_{QDH})\;.
\ee
The total energy $E$ to be minimized is then:
\be{total}
E = combine (E_{hard}, E_{medium}, E_{soft})\;.
\ee

\subsection{The generic Potts Algorithm}

The generic MFA algorithm appears in Figure~\ref{MFA-algorithm}.
Let ($N$) denotes the total number of the utilized neurons. At high
${T}$, the mean field solutions will tend to be states near to the
(symmetrical) maximum entropy state: $\forall (i, a) v_{ia} = {1/N}$.
Conversely at low ${T}$, finding a mean field solution will be
equivalent to using a local optimization method on the internal energy
--- a procedure highly sensitive to the initial conditions and known
to be ineffective (see p. 149 of \cite{32}).

These characteristics are similar to those of simulated annealing,
which is no surprise since both it and the mean field method compute
thermal averages over Gibbs distributions of discrete states, the former
stochastically and the later through a deterministic approximation.
It is therefore natural to couple the mean field method with the 
concept of annealing from high to low temperatures. The resulting
MFA method can be a powerful and general technique
which has found many applications in combinatorial optimization.

\begin{figure}[htb]
\vspace{0.25cm}
\begin{center}
\begin{enumerate}
\item Choose a problem $\Longrightarrow$ (mapped into) $\{w_{ij}\}$.
\item Find the approximate phase transition temperature by linearizing
  equation (\ref{form-5.2}).
\item Add a self-coupling $\beta$-term if necessary. Neural-wise, this
corresponds to a feedback connection from a neuron to itself.
\item Initialize the neurons $\nu_{ia}$ with ${1 \over N} \pm$ {\em
    random values}; and set $T(0) = T_c$.
\item Until ($\Sigma \approx {1 \over N}\sum_{ia}\nu^{2}_{ia}\geq 0.99$) do:
\begin{itemize}
  \item At each $T(n)$, update all $\nu_{ia}$ via $\upsilon_{ia}$:
  \begin{itemize}
    \item  $\upsilon_{ia} = - {1 \over T(n)} \cdot {{\partial
          E[\nu]} \over {\partial \nu_{ia}}}.$
    \item  $\nu_{ia} = {e^{\upsilon_{ia}} \over {\sum_b e^{\upsilon_{ia}}}}. \;\;\
; (\equiv f_{ia})$
  \end{itemize}
  \item $T(n+1) = 0.9 \times T(n)$.
\end{itemize}
\item After $\Sigma \geq 0.99$ is reached perform a {\it greedy
 heuristics} if needed to account for possible imbalances or rule violations.
\end{enumerate}

\end{center}
\caption{The Generic Mean-Field  Annealing Algorithm}
\label{MFA-algorithm}
\vspace{0.25cm}
\end{figure}

In addition to the structure of the energy function, there are three
major interdependent issues which arise in completely specifying 
a mean field annealing algorithm for a timetabling problem:

\begin{itemize}
\item The values of the coefficients of terms in the energy function.
\item The types of dynamics used to find solutions of the mean field
 equations at each ${T}$.
\item The annealing schedule details, i.e. the initial temperature
 ${T}(0)$, the rules for deciding when to reduce ${T}$ and
 by how much, and the termination criteria.
\end{itemize}

A quantity introduced by \cite{26} is called {\it saturation},
$\Sigma$, defined as
\be{saturation}
\Sigma = {1 \over N} \sum_{ia} \nu^{2}_{ia} \;\; ,
\ee
to characterize the degree of clustering\footnote{Other measures could
also do this, e.g. an entropy measure $-\sum_{ia}\nu_{ia}ln(\nu_{ia})$,
but $\Sigma$ will be used here for consistency with the original
work.} of events in time and/or space; 
clearly $\Sigma_{min} = {1 \over N}$ and $\Sigma_{max} = 1$.

\subsection{Some implementation details}

The general idea is to initialize the system close to the symmetric
state at a temperature near ${T}_c$, then, allowing the system
to converge to a fixed point at each temperature, anneal (i.e. reduce
${T}$) until either the system reaches a stable high saturation state
or some other termination criterion is satisfied.
Specifically, the implementation used here has the following
ingredients\footnote{These follow the method used in
\cite{26} but with some additions/improvements done by the
authors -see algorithm above.}: (Martin Simmen \cite{42}
investigated in depth these and most of the other issues pertaining to the Potts
model, using the TSP as a test case.)

Here are 4 of the main issues involved:

{\bf Initial State}: 
Another possible initial value of $\nu^{n=0}_{ia}$ is $(1 + \xi_{ia})/N$,
with the $\xi_{ia}$ being random variables drawn uniformly from 
${[-\xi, \xi]}$. Clearly $\xi$ should be $\ll 1$.

{\bf Convergence}:
At each ${T}$ the $\nu_{ia}$ variables are updated as above or by
the following rule\footnote{Thanks to M. Simmen for pointing this out.
Also see \cite{26}}:
\be{convergence-rule}
\nu^{(n+1)}_{ia} = \nu^{(n)}_{ia} + \gamma(f_{ia} - \nu^{(n)}_{ia}),
\;\; 0 < \gamma \leq 1\; .
\ee
Setting $\gamma = 1$ recovers the dynamics of \cite{26},
and the $\gamma \rightarrow 0$ limit gives the first order Euler
integration method for the differential equation:
\be{diff-equation}
d\nu_{ia} / dt = f_{ia} - \nu_{ia}\;.
\ee
The fixed points (also referred to as the free energy minima)
of this equation constitute the mean field solutions.
The degree of convergence is monitored through some sort
of a ``tolerance'' quantity of the form 
$\sum_{ia}{(\nu_{ia} - f_{ia})}^2$.

{\bf Annealing Schedule}:
A simple exponential scheme is used, i.e. after reaching a fixed point
at ${T}$, ${T} \mapsto {T} \cdot {T}_r$, where 
${T}_r \in (0,1)$ controls the cooling rate. Roughly speaking,
${T}_r$ should be set to $0.9$.

{\bf Termination}:
One criterion is as outlined above in the algorithm. The other two
involve terminating runs which either fail to reach a fixed point
within a certain number of sweeps at some temperature, or which allow
${T}$ to fall below some ${T}_{min}$, (${T}_{min} \ll
{T}_c$), without ever satisfying the $\Sigma > \Sigma^{thresh}$
criterion. ${T}_{min}$ was set to $10^{-3}$.

\subsection{More Constraints for the Potts  Model}

Let $q(i)$ represents a section of class $q$, $O_{q(i)}$ indicates
that section $q(i)$ is open, and $C_{q(i)}$ is closed. Also, let
$b(x,t)$ be a balanced interval for some time period $(t)$, $R$ is the
overall set of students, and an arbitrary given subset is  $r \subseteq R$. 

Here, we need to introduce the $R$-neuron: $S^{(R)}_{r(i);(d(t)d(x)}$ where
$k= 1,2,\ldots$.
Each student $r(i) \in R$ should have a balanced schedule 
over the two-group of days $(MWF)$ and
$(TuTh)$. We will also use $r_s$ to denote a single student.

For $X$ and $T$ neurons, let $x(j)$ and $t(j)$ where $j= 1,2,\ldots, n$
stand for places and time slots (or events), respectively. 
Using such notation, an example of the Potts normalization will be:
$S^{(X)}_{p(i)q(i);x(j)}$ and $S^{(T)}_{p(i)q(i);t(j)}$ respectively, where
$i=1,2,\ldots,n$ is a section index.

{\it The athletes case}: Let $\vec{r} \subseteq R$ denotes the group
of athletes  that are dealt with in this study. 
Let ($r^{*}$,$b$) be playing/training events scheduled during
the weekdays at (or during) times when academic (attending lectures)
events ($r^{+}$,$a$)  are also scheduled, or may possibly occur. Let
($d(x)$, $d(t)$) a place/time during the week at which any of these
events can occur. 

In addition to the above stated Potts neurons, we also will be using :

\begin{enumerate}
\item $S^{(\vec{r},X)}_{r^{*}b;d(x)}$ , $S^{(\vec{r},T)}_{r^{*}b;d(t)}$ :
It takes a value of ($1$) if event ($r^*,b$) takes place
in the place and time slot, respectively, of a weekday $d$, and zero otherwise.

\item $S^{(\vec{r},X)}_{r^{+}a;d(x)}$ , $S^{(\vec{r},T)}_{r^{+}a;d(t)}$ :
It takes a value of ($1$) if event ($r^+,a$) takes place
in the place and time slot, respectively, of a weekday $d$, and zero otherwise.

\item $S^{(\vec{r},X)}_{jk;d(x)}$ : 
It takes a value of ($1$) if event ($k$) of cluster ($j$) that is
associated with the $x$-neuron of $\vec{r}$ occurs in place $d(x)$,
and zero otherwise. 

\item $S^{(\vec{r},X)}_{\hat{j}\hat{k};d(x)}$ : 
It takes a value of ($1$) if event ($\hat{k}$) of cluster ($\hat{j}$)
that is associated with the $x$-neuron of $\vec{r}$ occurs in place
$d(x)$, and zero otherwise. 

\item $\tilde{S}^{(\vec{r},T)}_{j;d(t)}$ : 
It takes a value of ($1$) if cluster ($j$) that is associated with
the $t$-neuron of $\vec{r}$ occurs in time $d(t)$, and zero otherwise.

\item $\tilde{S}^{(\vec{r},T)}_{\hat{j};d(t)}$ : 
It takes a value of ($1$) if cluster ($\hat{j}$) that is associated with
the $t$-neuron of $\vec{r}$ occurs in time $d(t)$, and zero otherwise.
\end{enumerate}

In the above list, the last 2 neurons act as what Peterson calls {\it
effective} $time$-neurons. Not only will they handle periods of single
hours but also variants of longer and shorter intervals. These neurons
can be defined as follows:

Set ($\epsilon$) and ($c$) to some desired values depending on the time
period of the set of classes, now, an {\it effective}
$t$-neuron is:

\be{time-int.1}
\tilde{S}^{(\vec{r},T)}_{j;d(t)} = 
\sum^{g_j - c}_{k=0} S^{(\vec{r},T)}_{j;d(t-k+\epsilon)}
\ee

$\tilde{S}^{(\vec{r},T)}_{\hat{j};d(t)}$ can be defined in a similar manner.
If ($c = 1$) and ($g_i = 2$), then we are dealing with double-hour intervals
so we set ($\epsilon = 0$). If ($c = 1$) and ($g_i = 1$) then
($\epsilon = 0$) and

\be{time-int.2}
\tilde{S}^{(\vec{r},T)}_{j;d(t)} = S^{(\vec{r},T)}_{j;d(t)} .
\ee

Also, let $S^{(R)}_{q(i)r(i);b(x,t)}$ for balancing the enrollment in a multi-section
courses. Section $q(i)$ of students $r(i)$ are balanced over space-time $b(x,t)$.

To enroll student $r(i)$ in a particular section $q(i)$ given the tag $EC$ to
indicate the status of the student, 
we have $S^{(X[T])}_{q(i)r(i);(d(x)[d(t)])}|EC\in \{0,1\}$  is set to
$0$ if $EC$ is $0$, otherwise, it is set to $1$ and the action
$q(i)r(i)$ occurs during place $d(x)$ and/or time $d(t)$. 

Let $lin(c_i,c_j)$ representing the concept that course $c_i$ is linked to some
other course $c_j$ ; in this case, $c_i$ is a pre-requisite to $c_j$. 
The degree of ``linkage'' varies between $0$ and $1$. 
It takes a value of $1$ if under all circumstances, course $c_i$ must
be taken before $c_j$. To deal with this link concept, we assumed
above that we have an attribute $s(j,k)$ to represent a subject. For a
strong linkage, we will assume that we're dealing with subjects that
belong to the same cluster ($j$) (e.g. physics courses), so the link
between various events of ($j$) is: $lin(s(j,k), s(j,\hat{k}))$. 
Other degrees of linkage can involve subjects or courses from
different groups or clusters (e.g. physics and math courses, or 
engineering and math courses), and this connection is in the 
form: $lin(s(j,k), s(\hat{j},\hat{k}))$.

Let $\hat{x}(i)$ be a single location (or a seat) in class
$x(i)$ such that $i=1,2,\ldots, num$ ; where ($num$) is the 
number of class-rooms.

\subsubsection{Day-Class Balancing}
Each student of any particular section (i.e. $r(k)$) should have a
balanced schedule between ({\it MWF}) and ({\it TuTh}). The balancing
operation is accomplished by spreading as much as possible the courses
of the student over each of ({\it MWF}) and ({\it TuTh}).

Let $D$ be the set $\{\{MWF\}, \{TuTh\}\} = \{ D_1, D_2 \}$. Assume
the existence of the subject attribute $s(j,k)$. We introduce the
energy function $E_{QSD}$ that spreads the classes ($Q$) in a
particular field or subject ($S$) over the  set $\{D_1, D_2\}$.

The Potts condition of balancing $S^{({\bf R})}_{q(m)r(m);b(x,t)}$
(take a value of $1$ for balanced and $0$, otherwise) is accomplished
through one (or all) of these energies:

\be{balance.1}
E^{r}_{QSD} =
{1 \over 2}\sum_{j \neq \hat{j}}\sum_{k \hat{k}} 
\delta_{q(j,k),q(\hat{j},\hat{k})}
\delta_{s(j,k),s(\hat{j},\hat{k})}
\sum_{d_1,d_2} S^{(D)}_{j,d_1} S^{(D)}_{\hat{j},d_2}
\ee

\be{balance.2}
E^{r}_{QSD_1} =
{1 \over 2}\sum_{j_1 \neq \hat{j_1}}\sum_{k_1 \hat{k_1}} 
\delta_{q(j_1,k_1),q(\hat{j_1},\hat{k_1})}
\delta_{s(j_1,k_1),s(\hat{j_1},\hat{k_1})}\sum_{d_1}
\sum_{i \neq \hat{i}} S^{(D)}_{i,d_1} 
S^{(D)}_{\hat{i},d_1}
\ee

\be{balance.3}
E^{r}_{QSD_2} =
{1 \over 2}\sum_{j_2 \neq \hat{j_2}}\sum_{k_2 \hat{k_2}} 
\delta_{q(j_2,k_2),q(\hat{j_2},\hat{k_2})}
\delta_{s(j_2,k_2),s(\hat{j_2},\hat{k_2})}\sum_{d_2}
\sum_{i \neq \hat{i}} S^{(D)}_{i,d_2} 
S^{(D)}_{\hat{i},d_2}
\ee

\subsubsection{The Case of Student Athletes}
The main constraint of this case very much analogous to the {\bf hard} constraint
stated somewhere above : {\it different events should not occupy the same space-time slot}.

\be{athletic.1}
E^{\vec{r}}_{d(XT)} =
{1 \over 2}\sum_d\sum_{xt}\sum_{\dot{r}b}\sum_{\ddot{r}a} 
S^{(\vec{r},X[T])}_{\dot{r}b;d(x)[d(t)]} 
S^{(\vec{r},X[T])}_{\ddot{r}a;d(x)[d(t)]}
\ee

Where $(\dot{r}b)$ and $(\ddot{r}a)$ are some two different events. 
After some further refinement, we obtain:

\be{athletic.2}
E^{\vec{r}}_{d(XT)} =
{1 \over 2}\sum_{j \neq \hat{j}}\sum_{k \hat{k}}\sum_{d(x)} S^{(\vec{r},X)}_{jk;d(x)}
S^{(\vec{r},X)}_{\hat{j}\hat{k};d(x)} \sum_{d(t)} \tilde{S}^{(\vec{r},T)}_{j;d(t)}
\tilde{S}^{(\vec{r},T)}_{\hat{j};d(t)} + 
{1 \over 2} \sum_j g_j \sum_{k \neq \hat{k}}\sum_{d(t)}
S^{(\vec{r},X)}_{jk;x} S^{(\vec{r},X)}_{j\hat{k};x}
\ee

Since it is always the case that athletic training facilities 
are kept separate from the academic class-rooms and lecture halls, 
a more refined energy function is obtained:

\be{athletic.3}
E^{\vec{r}}_{d(T)} =
{1 \over 2}\sum_{j \neq \hat{j}}\sum_{d(t)} \tilde{S}^{(\vec{r},T)}_{j;d(t)}
\tilde{S}^{(\vec{r},T)}_{\hat{j};d(t)}
\ee

\subsubsection{Eligibility for Enrollment}
Before going through the process of balancing, each student 
($r(i) = r_s$) eligibility for enrolling in a given class need to be
determined.
Let the student's eligibility criteria be denoted by ($EC$), using 
the notation stated above, we could derive the
following energy function ---note: if $EC$ is set to ``No'' then it takes
a $0$ value in the energy function, otherwise, it is a $1$ :

\be{eligibility-criteria.1}
E_{QR} = 
{1 \over 2}\sum_{qr}\sum_{jk}\sum_{i}
\delta_{q,q(j,k)}\delta_{r,r(i)} S^{X}_{ij;x(i)}|({\it EC\in\{0,1\}})
\ee

\be{eligibility-criteria.2}
\sum_{i} O_{q(i)}=1
\Rightarrow \text{Test:}\left\{ EC(R) + \text{satisfy-constraints}(Hard,Soft)\right\}
\ee

\subsubsection{Course Pre-requisite and Subsequent Assignment}

Utilizing the link concept (${\bf lin(s,s')}$) between two different
subjects, a general energy assignment can be put together as follow: 

Let ($i$) stands for the situation at hand (also can be viewed as 
a space-time event).

\be{pre-requisite.1}
E_{link(1)} =
 -\gamma \sum_{is}{\left[\sum_{jk} \delta_{s,s(j,k)}
\tilde{S}^{(X[T])}_{ln(s,s');i}\right]}^{2}
\ee

\be{pre-requisite.2}
E_{link(2)} =
 -\gamma \sum_{s}\sum_{j \neq j'}\sum_{k \neq k'} \delta_{s(j,k),s(j',k')}
\sum_{i} \tilde{S}^{(X[T])}_{ln(s(j,k),s(j',k')); i}
\tilde{S}^{(X[T])}_{ln(s(j,k),s(j,k')); i}
\ee

\be{pre-requisite.3}
E_{link(3)} =
 -\gamma \sum_{j \neq j'}\sum_{k \neq k'} \delta_{s(j,k),s(j',k')}
\tilde{S}^{(X[T])}_{ln(s(j,k),s(j',k')); i}
\tilde{S}^{(X[T])}_{ln(s(j,k),s(j,k')); i}
\ee

If ($\gamma > 0$) then function minimization is the target here,
otherwise maximization will take place.

\subsubsection{Enroll in all parts  of a course}
There are a number course of different types that not only
have lecture periods but also associated sections that consist
of either labs or recitations or both. Signing up for such a course
entails enrolling in all associated sections of it. 

Let us assume the existence of some subject attribute $s(j,k)$
consisting of 2 components: $\dot{s}$ denoting lectures, and
$\ddot{s}$ stands for labs. Now the energy equation depicting
this assignment can be outlined as follow:

\be{enroll-in-all-parts.1}
E_{enroll(1)} =
{1 \over 2}\sum_{s}\sum_{jk}\sum_{i} \delta_{s,s(j,k)}\delta_{r,r(i)} 
\sum_{t} S^{(T)}_{jk;t}
\ee

\be{enroll-in-all-parts.2}
E_{enroll(2)} =
{1 \over 2}\sum_{si} \sum_{t\neq t'}\left[\sum_{jk} 
\delta_{s,\dot{s}(j,k)} \delta_{r,r(i)}
S^{(T)}_{jk;t}\right]\left[\sum_{jk'} \delta_{s,\ddot{s}(j,k')} \delta_{r,r(i)}
S^{(T)}_{jk';t'}\right]
\ee

\be{enroll-in-all-parts.3}
E_{enroll(3)} =
{1 \over 2}\sum_{k \neq k'}\delta_{s,\dot{s}(j,k)}
\delta_{s,\ddot{s}(j,k')}
\sum_{t \neq t'} S^{(T)}_{jk;t} S^{(T)}_{jk';t'}
\ee

\subsubsection{Lecture hall capacity}
In order to satisfy this constraint, the following inequality
need to hold.
\be{part-1}
\forall x \sum_{i}||x(i)| - \sum_{j}\sum_{k} S^{(X)}_{jk;x(i)}| \geq 0
\; .
\ee

\subsubsection{Student status and priority mapping}
\be{status-and-priority-mapping}
E_{assign} =
{1 \over 2}\sum_{m, \hat{x}, x}\sum_{i\neq i^{'}}
S^{(X)}_{\hat{x};x} 
\left( S^{(R)}_{r(i);(\hat{x}(m) \in x(m))}|({\it M = Req})
+
S^{(R)}_{r(i^{'});(\hat{x}(m) \in x(m))}|({\it M = \chi})\right)
\ee

\subsection{On the Energy Equation of the Potts Model}
The energy function of the Potts spin model is of the form:
\be{enrg}
E = E_{A} + E_{B} + E_{C} + \ldots
\ee

In \cite{26}, in which the TSP is used as study case,
parameters $\alpha$ and $\beta$ are used to control the strength of
E's restraints. The same parameters can also be used when dealing with
scheduling. 

\bo{ob-1}
\ \ For timetabling, $\beta$ is split into two parts ($\beta^{(T)}$
and $\beta^{(X)}$), one with respect to the time-based neurons and the
other with respect to the space-based neurons. $E_{\beta}$ becomes as
follows: 
\be{diagonal-energy}
- {\beta^{(T)} \over 2}\sum_j\sum_t {[S_{j;t}^{(T)}]}^2 - 
{\beta^{(X)} \over 2}\sum_{jk}\sum_{x} {[S_{jk;x}^{(X)}]}^2 
\ee
Here, the appearance of $\beta^{(T)}$ and $\beta^{(X)}$ is only in the
diagonal terms. They play a crucial role in monitoring the dynamics
(phase transitions, etc.). These diagonal terms correspond to self
coupling interactions for the neurons of the network. In our
simulations, each is set to $-0.1$.
\eo

\section{Some Performance Analysis of MFA}

Given $M_x$ = average number of rooms for each category of classes,
$N_x$ = the overall total number of rooms, and $N_t$ = total number of
slots , putting the figures of the students aside, with a computational
complexity of O($M_x N_x N_t$) per sweep, software implementation of
the MFA are unlikely to be competitive with conventional serial
algorithms such as 3-opt and Lin-Kernighan's for problems such as the
TSP, which empirically have an overall O($N^2$) running time
\cite{20}. Ideally, one would like to study the likely performance of
the method through numerical simulations using a small step-size
$\gamma$ to approximate continuous time dynamics. However, this is
computationally very expensive: not only do the simulation times grow
as $1 \over \gamma$, but also small values of ($\gamma$) accentuate
the difference between the expected value $T_c^{exp}$ of the critical
temperature $T_c$ and the theoretically computed value
$T_c^{theo}$. That is, low $\gamma$ runs tend to depress $T_c^{exp}$
artificially and thereby degrade solution quality, and as we observed,
result in a number of constraint violations. Many of those violated
constraints are hard, but the majority are classified to be soft to
medium. This later point can be compensated for by more stringent
annealing parameters, but at the cost of longer simulation times. No
attempts in our part have been made to carry this out. Also, we
have observed that the overall quality of the solution degraded as 
the problem size increased. This degradation spans both soft/medium
as well as hard constraints, and hence, illegal timetables.

%% From various correspondences with Andrew Gee, we would like to draw
%% the following as an observation:
%%
%% \bo{MFA-problem}\ \
%% We observe that the NP-completeness of the problem (the timetabling,
%% that is) is manifested in the fact that the constraint polytope
%% (\cite{14} and see previous disscusion on hypercube vertices) is
%% non-integral. Furthermore, for a decent number of constraints, there 
%% will be many (probably in the thousands) time more non-integral
%% vertices than 0-1 vertices. [Note here, since there are many objectives,
%% we just have to find any 0-1 point within the polytope to get a valid
%% timetable]. So, chances are that as we use some annealing mechanism to
%% derive the state vector towards the boundary of the polytope , it's
%% going to settle into a non-integral vertex, which is no use to anyone.
%% (Actually, Peterson alludes to these states in his paper with Gislen,
%% claiming that they are easily decomposed into two good solutions!)
%% Basically, there is no getting around this. Thinking of it from a
%% polyhedral \cite{14} perspective, we're trying to find the one integral vertex
%% in several thousands non-integral ones, and we're just driving the
%% state vector in a random direction towards the polytope boundary:
%% well, chances that it's not going to work.
%% \eo
%%
One thing that remains unclear is to what extent very slow annealing
can bring MFA close to finding the optimal (not only legal) schedule.
Whereas simulated annealing, given certain conditions, guarantees
(theoretically) convergence to the optimal solutions for asymptotically slow
annealing, such results is not true for MFA algorithms, because MFA is
an approximation.

\section{Experimental Analysis}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Our computations were done with a number of goals in mind. The main
objective was to provide a schedule which satisfied all hard
constraints and minimized the cost of medium and soft constraints, 
using real-life data sets for a large university. 
Our study case involved real scheduling data covering three semesters
at Syracuse University. The size and type of the three-semester data
is shown in Table \ref{table.01}. Nine types of rooms were used:
auditoriums, classrooms, computer clusters, conference
rooms, seminar rooms, studios, laboratories, theaters, and
unspecified types. Staff and teaching assistants are considered part
of the set of professors. Third semester (summer) data was much smaller
than other semesters, however, there were additional space and time
constraints and fewer available rooms. Our data was quite
large in comparison to data used by other researchers. For example,
high school data used by Peterson and colleagues \cite{27, 28} 
consists of approximately 1000 students, 20 different possible 
majors, and an overall periodic school schedule (over weeks). In the
case of Abramson {\it et al.} \cite{6}, their data set was created randomly
and was relatively small, and they stated that problems involving more
than 300 tuples were very difficult to solve. 

Table \ref{table.01} lists all major components of the data and Tables
\ref{table.02} and {\ref{table.03} show the rates of sparseness of the
three-semester data which greatly depends on number of time slots and 
rooms used. The level of sparseness increases as the overall size of
the data, particularly, number of rooms, time slots, and number of 
classes. Also, after scheduling all the classes, there will be about
$N_{sp} = (N_x N_t - N_l )$ spare space-time slots, hence, the
sparseness of the problem is defined as as the ratio $N_{sp} / (N_x
N_t)$. As we can see that Table \ref{table.03} degree of sparseness is
less than that of Table \ref{table.02} and the reason is the presence
of the students' preferences. Including student preferences makes the
problem much harder, but these are viewed as medium constraints and
thus are not necessarily satisfied in a valid solution.
 
Table \ref{table.1} shows the percentage of classes that could be scheduled in 
accordance with the hard constraints. Other constraints such as medium
and soft were also taken into account. Since the annealing process
in MFA is ``deterministic'' the results shown were averaged over 10
runs each having different initial conditions.
As expected, due to a higher sparseness ratio of the third (summer)
semester, the method does quite well. Also, we have observed that 
the overall results are generally below those of SA and ES
\cite{48}. In addition, we have found in the implementation of
this method that the results were quite sensitive to the size of the
data as well the type of constraints involved. If we confine ourselves
to the set of hard constraints, the results are as good as or even
better than the other methods. However if we take into account the
medium and soft constraints, that is, the overall cost function, this
method does not perform as well. 

\section{Conclusion}
We have applied the method of Mean Field Annealing to the difficult
problem of academic scheduling for a large university. We found that  
MFA works well for small scheduling problems, but does not appear to
scale well to large problems with many complex constraints. 
It is more difficult to tune the parameters, and use a preprocessor to
take advantage of a good initial solution, than methods such as
simulated annealing. Because of the complexity and size of the Potts
neural encoding in MFA, there seems to be no clear way of preserving
the state of an initial  configuration provided by a preprocessor when
carrying out the Potts encoding.

\newpage
\begin{table}[htb]
\vspace{0.1cm}
$$
\begin{tabular}{|l||c|c|c||}
\hline
                  & First Semester & Second Semester & Third Semester \\ \hline\hline
Rooms             &     509        &      509        &     120        \\ \hline
Classes            &     3839       &      3590       &     687        \\ \hline
Professors        &     1190       &      1200       &     334        \\ \hline
Students          &     13653      &      13653      &     2600       \\ \hline
Buildings         &     43         &      43         &     11         \\ \hline
Schools/Colleges  &     20         &      21         &     17         \\ \hline
Departments/      &     143        &      141        &     108        \\
Course Prefix     &                &                 &    \\ \hline
Areas of Study  (majors) &    200      &     200      &     200    \\ \hline
\end{tabular}
$$
\caption{Size of a three-semester data.}
\label{table.01}
\vspace{0.25cm}
\end{table}

\begin{table}[htb]
\vspace{0.1cm}
$$
\begin{tabular}{|l||c||}
\hline
Academic               & sparseness \\ 
Time Period           &    ratio   \\
                               &            \\ \hline\hline 
First Semester       &    0.50    \\ \hline
Second Semester &    0.53    \\ \hline
Third Semester     &    0.62    \\ \hline
\end{tabular}
$$
\caption{A three-semester sparseness ratio each (\underline{students
were not included}) computed using $N_{sp}/(N_x N_t)$ is shown. The
smaller the ratio, the less sparse the problem is and for this problem
the ratio is based on $N_t$ = 15 time slots each with a duration of 55
minutes. The first slot starts at 8:00 am and the last slot ends at
9:45 pm.}
\label{table.02}
\vspace{0.25cm}
\end{table}

\begin{table}[htb]
\vspace{0.1cm}
$$
\begin{tabular}{|l||c||}
\hline
Academic               & sparseness \\ 
Time Period           &    ratio   \\
                               &           \\ \hline\hline 
First Semester       &    0.35  \\ \hline
Second Semester &    0.38    \\ \hline
Third Semester     &    0.50    \\ \hline
\end{tabular}
$$
\caption{A three-semester sparseness ratio (\underline{including
students}) computed by ${\cal R}(N_{sp}/(N_x N_t)$, ((majors  *  max classes
taken per students) / number of students)). The smaller the ratio, the
less sparse the problem is and for this problem the ratio is based on
$N_t$ = 15 time slots each with a duration of 55 minutes. The first
slot starts at 8:00 am and the last slot ends at 9:45 pm.}
\label{table.03}
\vspace{0.25cm}
\end{table}

\begin{table}[htb]
\vspace{0.25cm}
$$
\begin{tabular}{|l||c|c|c||}
\hline
Academic        &   Scheduled  & highest   & lowest    \\
Time Period     &   (average)  & Scheduled & Scheduled \\
                         &    \%      &    \%     &    \%     \\ \hline\hline
First Semester        &       65.60  &   71.00   &   61.00   \\ \hline
Second Semester  &       67.20  &   75.00   &  65.00             \\ \hline
Third Semester      &      88.00   &  95.00    &  82.00 \\ \hline
\end{tabular}
$$
\caption{A three-semester MFA results: percentage of scheduled
  classes, averaged over 10 runs of different initial conditions.}
\label{table.1}
\vspace{0.25cm}
\end{table}

\clearpage
%\newpage
\begin{thebibliography}{HSW}

\bibitem{2}
Abramson, D,
``Constructing school timetables using simulated annealing: sequential
and parallel algorithms,'' {\it Management Science}, 37(1), 98-113,
1991.

\bibitem{6}
Abramson, Dang, and Krishnamoorthy,
``An Empirical Study of Simulated Annealing Cooling Schedules'',
Griffith Univ. report, Nathan, Qld, Aus. (1994).

\bibitem{39}
{\it Applied Simulated Annealing}, edited by Rene Vidal, Lecture
Notes in Economics and Mathematical Systems no. 396, Springer-Verlag,
1993.

\bibitem{58}
Carter, M, and C. Tovey, ``When is the classroom assignment problem
hard?,'' Working Paper 89-03, Dept. of Industrial Engineering,
Univ. of Toronto, 1989. (was submitted to {\it Operations Research}).

\bibitem{37}
Colorni, A, M. Dorigo, and V. Maniezzo, 
``On the use of genetic algorithms to solve the timetable problem'', 
report 90-060, Dipartimento Di Elettronica E Informazione, 
Politecnico Di Milano, Milano, Italia.

\bibitem{3}
Dowsland, K,
Using Simulated Annealing for Efficient Allocation of Students to
Practical Classes. Stat and OR Group, EBMS, Univ College of Swansea,
UK, 1993.

\bibitem{21}
de Werra, D,
``An introduction to timetabling,'' {\it European Journal of
Operational Research}, 19, 151-162, 1985.

\bibitem{49}
Eiselt H. A, and G. Laporte, ``Combinatorial Optimization Problems
with Soft and Hard Requirements,'' {\it J. Operational Research Society},
vol. 38, No. 9, pp. 785-795, 1987.

\bibitem{48}
Elmohamed, S, P. Coddington, and G. Fox,
``A comparison of Annealing Techniques for Academic Course
Scheduling'', submitted to the 2nd international conference on the
practice and theory of automated timetabling.

\bibitem{25}
Even, S., H. Itai, and A. Shamir,
``On the complexity of the timetable and multicommodity flow
problems'', SIAM Jr. Computing 5 (76), 691-703.

\bibitem{12}
Gary, M. and D. Johnson,
{\it Computers and Intractability: A guide to the theory of
NP-completeness}. Freeman, San Francisco (1979).

%% \bibitem{40}
%% Gee, Andrew,
%% Problem Solving with Optimization Networks, PhD Thesis, The Univ. of
%% Cambridge, UK, 1993.

\bibitem{45}
Gee, A, private communication.

\bibitem{27}
Gisl$\acute{e}$n, L, B. S$\ddot{o}$derberg, and C. Peterson, 
``Teachers and Classes with Neural Nets'', 
{\it Inter. Jr. of Neural Systems} 1, 167 (1989).

\bibitem{28}
Gisl$\acute{e}$n, L, B.  S$\ddot{o}$derberg, and C. Peterson, 
``Scheduling High Schools with Neural Nets'', Lund 
University Preprint LU-TP-91-9, Lund, Sweden (1991).
(appeared in {\it Neural Computation} 1992)

\bibitem{5}
Gosselin, K, and M. Truchon, ``Allocation of Classrooms by Linear
Programming'', {\it Journal of Operational Research Society} 37, 561 (1986).

\bibitem{38}
Hertz, A. ``Finding a feasible course schedule using Tabu search'',
{\it Discrete Applied Mathematics} 35, 255-270, 1992.

\bibitem{1}
Hertz,J., A. Krogh, and R. Palmer,
{\it Introduction to the theory of neural computation}, 
Addison-Wesley, Redwood City, 1991.

\bibitem{32}
Hopfield, J, and D. W. Tank, 
``Neural Computation of Decisions in Optimization Problems'',
{\it Biological Cybernetics} 52, 141 (1985).

\bibitem{33}
Hopfield, J. ``Neurons with graded response have collective
computational properties like those of two-state neurons''.
{\it Proceedings of the National Academy of Science } USA 81, 3088-3092
(1984).

\bibitem{20}
Johnson, D.
``Local Optimization and the TSP,'' Springer {\it Lecture Notes in CS} 443,
446-461, Springer-Verlag, 1990.

%% \bibitem{19}
%% Johnson, D,  C. Aragon, L. McGeoch, and Schevon,
%% ``Optimization by Simulated Annealing: an Experimental Evaluation, Part
%% I (Graph Partitioning).'' {\it Operations Research}  37, 865-892 (1989).

% \bibitem{10}
% Ladkin, P.B., and Maddux, R.D.,
% ``On binary constraint problems'', {\it Journal of the ACM} 41 (1994),
% pp. 435-469.

%% \bibitem{7}
%% Lister, R. ``Annealing Networks and Fractal Landscapes'', 
%% {\it Proc. IEEE Int. Conf. on Neural Nets}, March 1993, Vol. I, pp 257-262.

% \bibitem{11}
% Mackworth, A. K.,
% ``Consistency in networks of relations'', {\it Artificial
% Intelligence} 8 (1977) pp. 99-118. 

\bibitem{8}
Mouritsen, O. G,
Computer Studies of Phase Transitions and Critical Phenomena.
Springer-Verlag, Berlin, 1984.

% \bibitem{9}
% Montanari, U,
% ``Networks of constraints: fundamental properties and applications to
% picture processing'', {\it Information Sciences} 7 (1974), pp. 95-132.

%% \bibitem{46}
%% Peterson, Carsten, private communications.

\bibitem{29}
Peterson, C,
``Parallel Distributed Approaches to Combinatorial Optimization
Problems - Benchmark Studies on TSP'', {\it Neural Computation} 2, 261 (1990).

\bibitem{26}
Peterson, C, B. S$\ddot{o}$derberg,
``A New Method for Mapping Optimization Problems onto Neural Nets'',
{\it International Jr. of Neural Systems} 1, 3 (1989).

\bibitem{30}
Peterson, C, B. S$\ddot{o}$derberg, 
``Artificial Neural Nets and Combinatorial Optimization Problems'',
To appear in {\it Local Search in Combinatorial Optimization}, edited
by Aarts and Lenstra , Wily and Sons 1997.

\bibitem{50}
Schaerf, A, ``A survey of automated timetabling,'' Department of
Software Technology, Report CS-R9567, CWI, Amsterdam, The Netherlands.

\bibitem{14}
Shrijver, A,
{\it Theory of Linear and Integer Programming}, John Wiley and Sons,
Chichester, UK, 1986.

%% \bibitem{42}
%% Simmen, M,  ``Network Optimization'', PhD Thesis, Univ. of Edinburgh, 1993.

\bibitem{34}
Simic, P., ``Statistical mechanics as the underlying theory of
'elastic' and 'neural' optimizations'', {\it Network} 1, 89-103 (1990).

\bibitem{35}
Simic, P.,  ``Constrained Nets for Graph Matching and Other Quadratic
Assignment Problems'', {\it Neural Computation}, v.3, n.2, 268-281 (1991).

\bibitem{43}
P. Simic, private communication.

\bibitem{42}
Simmen, M, private communication.

%% \bibitem{4}
%% Thompson, J., and K. Dowsland, ``Variants of Simulated Annealing for
%% the Examination Timetabling Problem'', Working Paper, Statistics and
%% OR Group, European Business Management School, Swansea SA2 8PP, UK. (1994)

%% \bibitem{17}
%% Van Laarhoven, and Aarts,
%% {\it Simulated Annealing: Theory and Applications}. D. Reidel,
%% Dordrecht (1987).

\bibitem{18}
Van den Bout, D. E, and T. K. Miller III. Graph partationing using
annealed neural networks. {\it IEEE Transactions on Neural Networks},
1(2):192-203, June 1990.

%% \bibitem{41}
%% White, S. R. ``Concepts of scale in simulated annealing''. In
%% {\it Proceedings of the IEEE Int. Conference on Circuit Design}, pp 646-651,
%% 1984. 

\bibitem{36}
Wu, F, Y. ``The Potts model'', {\it Reviews of Modern Physics} 54,
235-268 (1982).
\end{thebibliography}

%% \appendix
%% \newpage
%% \renewcommand{\thesection}{Appendix A}
%% \renewcommand{\theequation}{\mbox{A}\arabic{equation}}
%% \pagestyle{myheadings}
%% \markright{Listing - Algorithm + Data Structures}
%% \setcounter{equation}{0}
%% \section{}
%% \appendix
%% \newpage
%% \renewcommand{\thesection}{Appendix B}
%% \renewcommand{\theequation}{\mbox{A}\arabic{equation}}
%% \pagestyle{myheadings}
%% \markright{Listing - the source code}
%% \setcounter{equation}{0}
%% \section{}

\end{document}

