%% 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 DRAFT-Report.97}
\end{flushright}

\vspace{1cm}
\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, Simulated 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
\begin{flushright}
{\sl DRAFT-Report.97}
\end{flushright}
\vspace{.5cm}
\begin{center} 
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\hrule
\vspace{.05cm}
\hrule
\vspace{.25cm}
{\bf 
Mean Field Annealing, Simulated 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}

\vspace{1.5cm}
\centerline{ABSTRACT}

In this study we tackle a complex problem of academic course 
scheduling (or timetabling) at the university level. First, we utilize
feed-back artificial neural system models along with the mean field theory 
approximation. We also extend the Potts normalization conditions and
equality constraints to accommodate our study case. The second
approach used is that of simulated annealing using two different annealing
schedules. A rule-based system is used as a pre-processor in both approaches.

\vfill
\noindent{\bf Key words : feed-back neural nets, MFT, 
tanh functions, Potts systems, constraints, events, stochastic
updating, annealing, simulated annealing, cost functions, temperature,
annealing schedule, phase transition, rule-based system.}

\vfill
\centerline{\today}
\vfill

\newpage
\pagestyle{empty}
\tableofcontents

\newpage

\pagestyle{plain}
\pagestyle{myheadings}
\markright{\sl DRAFT-Report.97}
\section{Introduction}
\setcounter{page}{1}

The prime objective of this study is to derive an approximate solution
to a high-dimensional and non-geometric combinatorial optimization (CO)
problem, which is academic scheduling at the university level. 
 
From the various categories of artificial neural systems (ANS)
for tackling CO problems, we choose one most suitable for
the problem at hand; which is the family of feed-back ANS based on
multi-state neurons with mean field theory (MFT) equations controlling
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. 

The second approach is that of simulated annealing. Three annealing
schedules were used: the standard geometric cooling schedule, the
 re-heating as a function of cost, and the adaptive annealing.

\section{Background and Motivation}

The problem of academic scheduling or timetabling has been tackled 
by many researchers \cite{47}, mostly in the field of
operations research. This is a very difficult problem due to its
nature as a multi-objective non-Euclidean combinatorial optimization
problem for which no good solution has previously been found.

Lately, there have been a number of {\it non-conventional} (or novel)
methods utilized to tackle this problem from different angles. Chief 
among them are neural networks \cite{32}, simulated annealing
\cite{2, 39}, genetic algorithms \cite{37}, 
tabu search \cite{38}, and various combinations of two or more of these
techniques. In this paper, we shall confine ourselves to simulated
annealing using the schedule of re-heating as a function of cost, and
mean field annealing. A rule-based system is used as a pre-processor
for each of these method.

Much of the work in the neural nets field 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. 

%% To achieve this, one needs a representation of the problem
%% in which the feasible solutions lie at ``hypercube'' vertices (more on
%% this later), and a quadratic energy function for which the minima
%% correspond to solutions and the depth of each minimum reflects the
%% solution quality. 

In practice however, considerable difficulties have been encountered
carrying out this program for the TSP, leave alone non-Euclidean CO
problems such as scheduling. 

\section{The Timetabling Problem - description and solution methods}
One of the problems that does not have a well-defined objective
function to optimize and of the class of constraint satisfaction
problems is the $NP$-complete  {\it teachers and classes} timetabling
problem studied in \cite{2, 3, 27}.

Briefly, the problem is defined as follows:

\bd{general-TTS}\ \
For a certain school with $N_p$ teachers, $N_q$ classes, $N_x$
classrooms or lecture halls, and $N_s$ students, it is required to
schedule $N_l$ lessons (defined as teacher-class pairs) within a time
limit of $N_t$ lesson periods (or time slots). 

A legal schedule needs to be found such that no teacher, class,
or student is in more than one place at a time, and no lecture hall is
expected to accommodate more than one lesson at a time or more
students than its capacity. 
\ed
From this point on, we will refer to academic scheduling or
timetabling as $TTS$ (short for a timetabling system).


\subsection{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}.

\bd{CSP}\ \
{\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}
\ed

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.

\subsection{Constraints of TTS}

In a nutshell, timetabling is the assignment of time slots to a set of
events, {\bf subject to} constraints on these assignments. The medium
and soft constraints have an associated cost and if they are not
satisfied, the goal is to minimize this cost. In regard to the hard
constraints, one must derive their associated cost (or penalty) to
zero.

These constraints are divided into three types:
\begin{enumerate}
\item Hard - These constraints that must be satisfied.
\item Medium - These constraints should be satisfied, if possible, but
are not required for a valid solution. 
\item Soft - These constraints have a lower priority but always the
preference is to minimize them.
\end{enumerate}

\subsubsection{Hard Constraints}

The most common and important hard constraint is what is sometimes
referred to as an {\it edge} constraint between two events. This
reference is due to the similarity of any {\bf simple} timetabling
problem with {\it vertex coloring} in graph theory; which simply
states that any given pair of events ${\bf e_1}$ and ${\bf e_2}$ must
not overlap in time. ********[more on this, shortly]************

Those events that must not overlap in time include:

\begin{itemize}
\item two classes taught by the same professor,
\item a class and a recitation or a lab of the same class,
\item classes with students in common.
\end{itemize}

Other general hard constraints include:

\begin{itemize}
\item A class can not be assigned to a particular room unless {\it
size of the room} $\geq$ {\it the class enrollment}. The smaller the
difference between these two sizes, the better the assignment.
\item A limit on how many classes taught by each instructor is enforced.
\end{itemize}

Other constraints are those directly related to students such as:
\begin{itemize}
\item To sign up for a class or section, it must be open.
\item Eligibility criteria for the section must be met, etc.
\item Do not enroll athletes in courses that conflict with their sport
practice time, etc.
\end{itemize}
These constraints have to be respected at all costs.

\subsubsection{Medium Constraints}

These kinds of constraints fall into, so to speak, the shaded region
between the {\it hard} and {\it soft} constraints in the sense that
they do not deal with time conflicts but do have a high penalty
attached to them. Also, this penalty is not as high as that associated
with the hard constraints. In the final (legal) schedule the penalty
of these constraints is minimized and it is quite preferable to reduce
it to $0$.

In short, these kind of constraints should be satisfied but are
nevertheless relaxable.

{\bf\underline{Examples of medium constraints}}
\begin{itemize}
\item For each student, balance {\it MWF} and {\it TuTh} sections.
\item Balance or spread out the lectures over the week.
\item Classes may request contiguous time slots.
\item Classes may be grouped with other classes.
\item The overall schedule should have as few ``holes'' as possible;
that is lessons should be glued together, etc.
\end{itemize}

\subsubsection{Soft Constraints}

Like medium constraints, they also do not deal with time conflicts but
unlike the medium constraints, have a low penalty (or cost) associated
with them. Furthermore, this cost is lower than that of the medium
constraints and much lower than that of the hard ones. In addition, 
in the final schedule, this cost is minimized but does not often
approach zero.

{\bf\underline{Examples of soft constraints}}
\begin{itemize}
\item Balance enrollment in multi-section courses.
\item Lunch and other break times may be specified.
\item Teachers may request periods in which their classes are not taught.
\end{itemize}

\subsection{The overall cost function}

The {\it cost function} measures the infeasibility of the current
solution and may involve the weighted sum of penalties associated with
different types of feasibility violation.

We consider the following as part of the function:
\begin{itemize}

\item Distance cost: proportional to the distance between
where the class is assigned and the building which is housing its home
department or school.

\item Time period cost consisting of :
\begin{enumerate}
\item class cost, 
\item professor cost, 
\item cost of class spreading over the week and within each of the
weekdays, 
\item cost of class glueing and filling in any holes of the timetable,
\item room size cost, and 
\item students cost which is the sum of various costs associated
with each student, such as eligibility criteria and number of classes
taken by the student.
\end{enumerate}

\item The cost of professor is computed from:
\begin{enumerate}
\item cost incurred from  assigning more classes than the
maximum number taught by that professor, plus 
\item cost from not having a long enough time interval between courses
taught by the same professor, and
\item cost of assigning to the same professor two or more classes that
overlap in time. 
\end{enumerate}

\item The cost of room size is computed in a non-linear fashion
depending on the magnitude of the difference between the class
enrollment and the room size.
\end{itemize}

{\bf Note:} {\em Using summation is not strictly necessary, provided that
the components can be combined in some fashion to reflect the quality
of the solution.}
  

\subsection{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-complete problem, see Gary and Johnson
\cite{12}.

The problem 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.

Using the LP representation, for a set as small as 20 teachers,
30 time-intervals, and 10 classes, the problem would be represented by
6000 variables. This makes the problem intractable.

The overall complexity of timetabling is in the order of 
$2^{{(N_c N_x N_t)}^a}$, for some  $a \ge 1$.

For a much simpler problem, which is the {\it class allocation} (no
timetabling part included), using linear programming
\cite{5}, the computational load ($n_c$) scales as follows:

\be{complexity}
n_c \propto N^3_c N^2_x N^3_t\;\;,
\ee

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


\section{Some Methods of Tackling Timetabling}

%%
%% From the work of Peterson \cite{Peterson-89.2} and Gee
%% \cite{Gee-93}, the problem is mapped using a 0-1 $N_p N_q N_x
%% N_t$-element\footnote{Briefly, the solution of this CO problem involves the
%% minimization of one or more energy (or cost) function(s) -other
%% problems, such as the TSP, involves the minimization of a single
%% quadratic objective function-  over a set of variables between $0$ and
%% $1$. The elements of $v$ range over the 0-1 interval.} solution vector
%% ($v$), which describes how teacher-class pairs (i.e. lessons) are to
%% be associated with space-time slots. 
%%

\begin{figure}[htb]
\centerline{\psfig{figure=polytopes-2.ps,width=6in,angle=-90}}
\caption{An integral polytope has all its vertices at integral points
(points whose coordinates are all integers). It follows that an
integral polytope within the unit hypercube has all its vertices at
$[0-1]$ points.}
\label{poly}
\end{figure}

The timetabling problem has no cost function, only constraints (hard,
medium and soft) each with their associated penalty (or cost); all we
have to do is find any integral vertex within the constraint polytope
$P$ (see Figure \ref{poly}). To elaborate a bit on this, using the 
work of  Schrijver \cite{14, 15}: we think of each constraint as a
hyperplane in the solution space. The set of constraints define a
bounded polyhedron (or polytope) within which the state vector must
remain if it is not to violate any constraints. The $0$ or $1$
solution points will lie at some of the vertices of this polytope, and
whatever annealing method is used, it must force convergence to one of
these points. The problem with annealing (of course, depending on the
schedule used) is that the best one can expect is convergence to a
vertex, not necessarily a 0 or 1 vertex. In the case of MFA, the
ability of the neural net finding a valid solution point reliably
therefore depends on the polytope having vertices of 0 or/and 1 values
(in short, 0-1 vertices), in which case it is called an integral
polytope. Now, polytopes like the one for the permutation matrix
mapping for the TSP are integral, so there is no problem finding valid
TSP solutions (though the tour length is another question
altogether). Class scheduling polytopes are a totally different story:
Gee \cite{40} and some of his colleagues at Cambridge University
have done some experiments on them and found that for timetabling
problems as small as five teachers, five classrooms, five classes, and
five time-slots, the proportion of 0-1 vertices is only about
$1\%$. This mean that on average one has to run a descent procedure
about 100 times to find a valid solution!
Some of our experiments with sets of teachers, classes, classrooms,
students and time-slots of various sizes confirm some of Gee's results
and conclusions. 

Timetabling problems can be characterized by their {\it sparseness}.
After all the required lessons have been scheduled, there will be $N_s
= (N_x N_t - N_l)$ spare space-time slots, hence, the sparseness of
the problem is defined as the ratio $N_s / (N_x N_t)$. The sparseness
decreases as  $N_s / (N_x N_t)$ goes to zero or as $N_s$ becomes
smaller and smaller. Furthermore, \cite{40} notes that for sparse
problems, the probability of an arbitrary descent procedure converging
to an integral vertex is near one, though this falls off rapidly for
less sparse problems. For dense problems, there is an additional
correlation involving the problem size: as the problem gets larger
(more classes, teachers, rooms, students, etc.), the proportion of
integral vertices gets smaller. From this we could deduce that a
one-shot descent procedure, operating on an arbitrary ``objective
function'', is not likely to be very successful with difficult
constraint satisfaction problems, most probably converging to a
non-integral vertex of $P$, for example. On the other hand, any
``partial'' or ``imperfect'' solution to a timetabling problem found
by methods such as simulated annealing (using some good adaptive
schedules) is likely to be of some use. Since simulated annealing does
not search over the interior of the unit hypercube, but only over its
vertices, the solution it finds will be a unique timetable, possibly
-as it is the case in our experiments- violating as few as one or two 
constraints. As part of the post-processing operation, we could at 
least remove some lessons to obtain a valid (or legal) timetable which
fulfills most of the teaching requirements. Furthermore, we could do
few more swaps (or moves) of classes around to arrive at a timetable
meeting all the original requirements.

\section{Mapping of Timetabling to a Graph Coloring Problem}

\subsection{Graph Coloring}

Given a graph $G=(V,E)$, the graph coloring problem is to find a
coloring of its vertices with as few colors as possible. The main
constraint is that two vertices linked by an edge should not receive
the same color (such coloring are called feasible coloring). One could
define the set $S$ of feasible solutions as the set of feasible
coloring which do not use more than a given number of colors; for
example, this upper bound may be the number of colors used in the
initial solution.

\subsection{The Mapping}

\begin{enumerate}
\item The courses are the vertices of a graph.
\item The time periods correspond to the colors.
\item Two vertices are linked by an edge if the corresponding courses
cannot be given at the same time.
\item Additional constraints have to be taken into account and mapped
onto the graph.
\end{enumerate}

\subsection{Interpretation of the Solution Space}

After setting up or mapping the problem we need to determine the
following: 
\begin{itemize}
\item 
{\underline{A feasible solution}} is defined as a schedule that
satisfies a subset ($C^{'}$) of constraints ($C$) : $C^{'} \subseteq C\;$;
provided that all hard constraints are part of $C^{'}\;$.
\item 
A neighboring solution is obtained by modifying the current solution
space (i.e. swapping two or more classes in time).
\item
Here are two guidelines in determining whether a given constraint to be
included in $C^{'}$ or not:
\begin{enumerate}
\item A class should never be assigned to a time period if such an
assignment induces only class schedules that do not satisfy all the
constraints (specifically the hard constraints).
\item A conflict between two or more classes should not be a
sufficient condition for hindering an assignment.
\item
If all courses last one time period, this is equivalent to changing
the color of a vertex, for example.
\end{enumerate}
\end{itemize}


\section{The Preprocessor Stage}

\subsection{The Problem}

Given data sets of courses and their days, hours, enrollments, and
instructors; rooms and their capacities, specification, and locations;
distances between buildings; priorities of each building; students and
their courses: the problem is to construct legal school timetable(s)
satisfying all the hard constraints and minimizing the medium and soft
constraints. 

Courses that are not assigned to rooms due to conflict in time or size
are left out of the final timetable for the user to handle manually. 

\subsection{The Method}

The technique used is a rule-based system with lots of heuristics and
conventional recursion to assist in carrying out course
assignments. Here are few examples of some of the implemented rules
dealing with {\bf space and hours}, {\bf type of room}, and {\bf
priority of room}:

\subsubsection{Space and Hours}

{\it IF [room(capacity) $>$ course(space-requested)] and [no time conflict in
   this room]
THEN
   assign the room to the course:}

[If more than a room having an identical availability time are present
in the database, then choose the one with the minimal magnitude of the
difference between room(capacity) and course(space-requested)].

{\bf\underline{Input:}}
CPS-615 requests space for 75 students and meeting hours of [MWF 10:40 - 11:35 am]. 

In the database we have, for example, the following rooms:

\begin{itemize}
\item room: 1-218 CST , capacity: 45 students, and pre-occupied hours: none.
\item room: 207 HL, capacity 120 students, and pre-occupied hours: [MW
09:45 - 11:10 am]. 
\item room: Kitt-Aud HBC, capacity: 127 students, and pre-occupied
hours: [TTH 10:40 - 12:00 am].
\end{itemize}

{\bf\underline{Output:}}
Room Kitt-Aud is assigned to CPS-615. \ \

\subsubsection{Type of Room}

{\it IF [room(capacity) $>$ course(space-requested)] and [no time conflict in
   assigning this room] and [type(room) matches type(requested)] 

THEN
   assign the room to the course:}
[If identical rooms in type are
present in the database, then choose the one with the minimal
magnitude of the difference between room(capacity) and
course(space-requested)].

{\bf\underline{Input:}}
BIO-615 requests space (not biology lab) for 35 students and meeting
hours in [MW 10:40 - 11:35]. 

In the database we have, for example, the following rooms:

\begin{itemize}
\item room: 1-218 CST, capacity: 45 students, type: classroom, and
pre-occupied hours: none.
\item room: 124-SIMS, capacity: 40 students, type: laboratory, and
pre-occupied hours: none.
\end{itemize}

{\bf\underline{Output:}}
Room 124-SIMS is assigned to BIO-615. \ \

\subsubsection{Room Priority}

Rooms are classified as:
\begin{enumerate}
\item General Purpose: used by all users.
\item Priority Rooms: are used based on certain priority of the class
or the department.
\item Special Purpose: are only assigned to specific users.
\end{enumerate}

\begin{itemize}

\item {\it IF [room(priority) == general] and [room(capacity) $>$
course(space-requested)] and [no time conflict in assigning this room] and
[type(room) matches type(requested)]  

THEN  assign any user to this room.}

\item {\it IF [room has a priority user] and [the user is already
assigned] and [room(capacity) $>$ course(space-requested)] and [no
time conflict in assigning this room] and [type(room) matches
type(requested)]

THEN  assign other users.}

\item {\it IF [room(priority) == special purpose] and [user belongs to
the special users of this room] and [room(capacity) $>$
course(space-requested)] and [no time conflict in assigning this room]
and [type(room) matches type(requested)]

THEN  assign the room to the user.}

\item {\it IF [room has a priority user] and [the user is already
assigned] and [room(capacity) $>$ course(space-requested)] and [no
time conflict in assigning this room] and [type(room) matches
type(requested)]

THEN  assign other users.}

\item {\it IF [room(priority) == special purpose] and [user belongs to
the special users of this room] and [room(capacity) $>$
course(space-requested)] and [no time conflict in assigning this room]
and [type(room) matches type(requested)]

THEN  assign the room to the user.}
\end{itemize}

\subsection{Components of the Rule-Based System}

This system consists of a number of heuristics implemented to take
care of various cases that arise when carrying out class
scheduling. We have developed this system specifically for the problem
of academic scheduling. The basic data structures or components of the
system are:

\begin{itemize}
\item {\bf Distance Matrix} of integer unit values between each
academic department and every other building under use for scheduling.

\item {\bf Course Structure} of each course scheduled anywhere in
campus. These structures are capable of linking with each other.
 
\item {\bf Room Structure} of each room (regardless of type) involved
in the scheduling process. Like courses, rooms are capable of linking.

\item {\bf Time Periods Structures} to keep track of which hour or
time slot is occupied and which is not.

\item {\bf Structure of department inclusion} many departments belong
to other larger departments or colleges, etc.
\end{itemize}

\subsection{The Overall Function of the System}

\begin{itemize}
\item It takes in data from external sources:
\begin{itemize}
\item Course data supplied by the university.
\item Building-room data also provided the registrar's office.
\item Department-to-building distance matrix.
\item Inclusion data, which a list of all colleges and schools and
their associated departments or sections. Each list is headed by a
college or a school prefix followed by its associated departments
prefixes (i.e. ECS CS MAE  etc.).
\end{itemize}
 Then builds an {\it internal database}. 
\item Using this internal database, the scheduling process takes place.
\end{itemize}

The scheduling process involves:

\begin{itemize}
\item Checking distance.
\item Checking room type, building, and hours occupied.
\item Check and compare hours for any possible conflict.
\item Check rooms for any possible space conflict.
\item Keep track of and update the list of hours already scheduled 
      (or occupied). 
\end{itemize}

\subsection{The Output of the Rule-Based System}

The {\it first output} of this system consists of classes, with their
associated professors and students, assigned to various rooms.

The {\it second output} is a list of classes which could not be
assigned due to time/space or other conflict(s).

These two outputs are fed to the main heuristic which can be either a
{\it simulated annealing} or a {\it mean field annealing} algorithm to
yield a final legal schedule.

\subsection{The Rationale of Using a Preprocessor}
\begin{itemize}
\item It is a good practice when using heuristics such as
simulated annealing, to start at a good point in the search
space or with a partial solution of the problem; regardless of the
difficulty of the problem.

\item A good starting point can only have a positive effect on the
overall solution, and that is improving its quality.

\item A good starting point or a partial solution will also have a
positive effect on the convergence of the system; faster convergence
to the final solution.
\end{itemize}

\newpage
\subsection{Mapping the Problem Components onto Time Periods}

\begin{figure}[htb]
\begin{center}
\begin{picture}(380,380)
\psfig{figure=tuples-periods.ps,height=380pt,width=380pt}
\end{picture}
\end{center}
\caption{This shows how the main problem components are mapped onto
the desired time periods.}
\label{components-mapping}
\end{figure}

\newpage
\section{Mean Field Annealing}

\vspace{0.1cm}
\begin{table}[htb]
$$
\begin{tabular}{||ccc||ccc||}
\hline
 & Physical System    && Optimization Problem   &&\\ \hline\hline\hline\hline
 &  System State      &&  Feasible Solution     &&\\ \hline
 &  Change of State   &&  Neighboring Solution &&\\ \hline
 &  Energy            &&  Cost                  &&\\ \hline
 &  Ground State      &&  Optimal Solution      &&\\ \hline
 &  Rapid Quenching   &&  Local Search          &&\\ \hline
 &  Careful Annealing &&  Simulated Annealing   &&\\ \hline
 &  Temperature       &&  Control Parameter     &&\\ \hline\hline
\end{tabular}
$$
\caption{The analogy between physical systems and the optimization problems.}
\label{table.0}
\end{table}
\vspace{0.5cm}

To start with, a summary of the Hopfield network is given as a prelude
to mean field annealing.

\subsection{The Hopfield Net - A summary of the continuous model}

\begin{itemize}
\item The network is built by connecting a large number of neurons to
each other to form a dense, parallel array.
\item In general, any $i^{th}$ neuron in the net is described by two
variables: 
\begin{itemize}
\item its current state (input) $u_i$,  and
\item its output $v_i$.
\end{itemize}

\item The output is related to the input by a monotonically increasing
transfer function: $v_i = g_{\lambda}(u_i)$.

\item The net's dynamics are governed by the first order differential
equation:
\be{FODE.1}
{{\bf du} \over {\bf dt}} = 
- \eta {\bf u} + {\bf W v} + {\bf I^b}\;.
\ee

\item $g_{\lambda}(u_i) = {1 \over {1 + e^{-2\lambda u_i}}}\;$, for
large $\lambda$, this function approximates a step function,
forcing $v_i$ to $[0,1]\;\Longrightarrow\;$ (the output of the net
lies on the unit hypercube).

\item $W_{ij}$ is the connection strength between neuron $i$ and neuron
$j$. 

\item $I_i^b$ is an input bias to neuron $i$.

\item Neurons may be updated individually at random (asynchronously)
or all together (synchronously).
\end{itemize}

\subsubsection{Operation of the Network}

\vspace{0.1cm}
\begin{figure}[htb]
\begin{center}
\begin{picture}(320,325)
\psfig{figure=hopfield_net.eps,height=320pt,width=325pt}
\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:

\begin{itemize}
\item The network's 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}
{{\bf du} \over {\bf dt}} = - \eta {\bf u} + {\bf W v} + {\bf I^b}\;.
\ee

\item The $I_{i}^b$ corresponds to an input bias to neuron $i$.
\item 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.
\end{itemize}

\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}
{\bf \dot{u}} = 
-\eta{\bf u}+{\bf W v}+{\bf 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 {\bf 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 {\bf 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

\begin{itemize}
\item $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.

\item 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$.
\end{itemize}

\subsubsection{Updating the States of the Net}
\begin{enumerate}

\item Updating of the states (or the output of the neurons) is carried
out either asynchronously or synchronously. 
\item From a simulation point of view, we found that synchronous updating
a bit more efficient.
\item 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.
\end{enumerate}

\subsubsection{Mapping the Problem onto the Network}
\begin{itemize}
\item A representation of the problem is needed in which the feasible
solutions lie at hypercube vertices.

\item Each hypercube vertex is regarded as a configuration of the
problem, with an associated energy (or cost) function.

\item The partition function is then formed by summing the Boltzmann
factors of the admissible configurations.

\item 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.

\item Hopfield and Tank example is the TSP, which they have mapped
onto the net using a permutation matrix representation of the TSP
tours.
\end{itemize}

\subsubsection{Shortcomings of the Hopfield Model}
\begin{itemize}
\item The net's partition function sums over a vast number of
configurations (e.g. for the TSP, all of the ${2^{N}}^2$
configurations are admissible). 
\item 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.
\end{itemize}


\subsection{Mean Field Annealing: The Potts Approach 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 and S$\ddot{o}$derberg
\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 (definition \ref{norm-dfn})
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, 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, saddlepoint
equations are derived (see next sections), the solutions of which pick
out the dominant states of the network at the current temperature
${\cal 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 \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
${\cal 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{Notation and Definitions (The First Part)}

In this section we define some key terms and concepts used 
in this paper and primarily adapted from the work 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$.

\bd{dfn-0}
\ \ 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
\ed

\bd{const-dfn.1}
\ \ There are two types of constraints to deal with in this paper:
\begin{enumerate}
\item Those involving professors, classrooms, subjects, and time
periods.
\item Those dealing with students scheduling.
\end{enumerate}
In addition, some of these constraints are usually classified to be
{\bf hard} and the rest are {\bf soft}.
\ed

\bd{const-dfn.2}
Here, we define a number of {\bf hard constraints} of the first
type (see definition \ref{const-dfn.1}) .
\begin{itemize}
\item An event $(p,q)$ should occupy precisely one space-time slot $(x,t)$.
\end{itemize}
This 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
\begin{itemize}
\item Different events $(p_1,q_1)$ and $(p_2,q_2)$ should not occupy the same
space-time slot $(x,t)$.
\end{itemize}
\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
\begin{itemize}
\item A teacher $p$ should have at most one class at a time.
\end{itemize}
\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
\begin{itemize}
\item A class $q$ should have at most one teacher at a time.
\end{itemize}
\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
\ed

 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 
 ${\cal T}$ to denote the ${\it 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 {\cal T}}\right)
 \ee
When ${\cal 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.  


\bd{don-1}
\ \ At a finite ${\cal 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 {\cal 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$ Pr(neuron(i) $\in$ state(a)).
Equations (\ref{form-5.1}, \ref{form-5.2}) are computed in an iterative manner.
\ed

\bd{norm-dfn}
\ \ 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.
\ed

\bd{dfn-2}
\ \ To carry out the minimization of the energy functions of definition (\ref{dfn-0})
we need to convert the neural variables $\{0,1\}$ of $S_i$ to their
corresponding {\it mean field variables} at ${\cal T}$. These
correspondences are derived (see \cite{30}) to be: 

\be{form-6.1}
V_{pq;xt} \approx <S_{pq;xt}>_{\cal 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 {\cal 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
\ed

The $MFT$ equations of $\upsilon$ and $\nu$ are
solved {\it iteratively} either synchronously or asynchronously, under
${\cal 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}$ :

\bd{dfn-3} \ \
\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

\ed

During annealing, as we go from a higher temperature to a lower one, 
a {\it phase transition} is passed at (${\cal T} = {\cal T}_c$), where 
${\cal 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 ${\cal T}_c$ set the scale or the range of the value of
${\cal 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 ${\cal T}$ = 0} .
\ee

Also ${\cal 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}.

\bd{dfn-4}\ \
{\em Relative Clamping}: 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 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.
\ed

\bd{dfn-3.9}
To deal with events that take place over two consecutive hours, what \cite{26}
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.
\ed

Utilizing definition (\ref{dfn-3.9}), constraints (\ref{form-3},
\ref{form-4}, \ref{form-5}) can be put into the following form:

\bd{dfn-3.9.1}

\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
\ed

\bd{dfn-6}\ \
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
\ed

\bd{dfn-3.9.2}\ \
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 {\it quasi}-classes. This is a common practice
with those classes where the students have many choices.
To carry out such an extension to quasi-classes 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
\ed

\bd{dfn-3.9.3}\ \
From definitions ($\ref{dfn-4}$, $\ref{dfn-3.9.2}$), 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$ .
\ed

\bd{dfn-3.9.4}\ \
Utilizing the concepts of the above definition ($\ref{dfn-3.9.2}$,
$\ref{dfn-3.9.3}$), the collision terms of definition
($\ref{dfn-3.9.1}$) 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).
\ed

The factorization idea of (dfn \ref{dfn-3.9.4}) 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 above equations
(\ref{form-9.1}, \ref{form-10.1}) .

Substituting $V$ for $S$ in (\ref{form-9.1}, \ref{form-10.1}), 
the {\it mean field} equations at a finite 
${\cal 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:

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

\be{form-4.1.4}
U^{(T)}_{j;x} =
- {1 \over {\cal 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}} - {\ref{form-4.1.6}})
can either be {\bf synchronous} or {\bf serial} updating.
\ed

\bd{dfn-4.2}
{\em Synchronous Updating}: \ \
A single sweep consists of first computing all the local fields, and then
updating all the neurons. That means do equations ({\ref{form-4.1.3}},{\ref{form-4.1.4}})
then perform equations ({\ref{form-4.1.5}},{\ref{form-4.1.6}}).
\ed

\bd{dfn-4.3}
{\em Serial Updating}: \ \
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}}).
\ed

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

\bd{dfn-5}\ \
\begin{itemize}
\item {\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.
\end{itemize}

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

\begin{itemize}
\item {\em Glueing}: The schedule should have as few ``holes'' as possible; that is
lessons should be glued together.
\end{itemize}
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}, 
\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
\ed

\subsection{The generic Potts Algorithm}
\subsubsection{Some general comments on the method}

Let ($N$) denotes the total number of the utilized neurons. At high
${\cal 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 ${\cal 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
{\it mean field annealing} method can be a powerful and general technique
which has found many applications in combinatorial optimization.

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 ${\cal T}$.
\item The annealing schedule details, i.e. the initial temperature
 ${\cal T}(0)$, the rules for deciding when to reduce ${\cal 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$.

\subsubsection{The Algorithm}
\begin{itemize}
\item Choose a problem $\Longrightarrow$ (mapped into) $\{w_{ij}\}$.
\item Find the approximate phase transition temperature by linearizing equation 
(\ref{form-5.2}), see also above section.
\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
${\cal T}(0) = {\cal T}_c$.
\item Until ($\Sigma \approx {1 \over N}\sum_{ia}\nu^{2}_{ia}\geq 0.99$) do:
\begin{itemize}
  \item At each ${\cal T}(n)$, update all $\nu_{ia}$ via $\upsilon_{ia}$:
  \begin{itemize}
    \item  $\upsilon_{ia} = - {1 \over {{\cal 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 ${\cal T}(n+1) = 0.9 \times {\cal 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{itemize}

\subsubsection{Some implementation details}
The general idea is to initialize the system close to the symmetric
state at a temperature near ${\cal T}_c$, then, allowing the system
to converge to a fixed point at each temperature, anneal (i.e. reduce
${\cal 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.)

\begin{itemize}
\item {\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$.

\item {\bf Convergence}\ \
At each ${\cal 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$.

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

\item {\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
${\cal T}$ to fall below some ${\cal T}_{min}$, (${\cal T}_{min} \ll
{\cal T}_c$), without ever satisfying the $\Sigma > \Sigma^{thresh}$
criterion. ${\cal T}_{min}$ was set to $10^{-3}$.

\end{itemize}

\section{Our study case using the Potts model}

\subsection{Notation and Definitions (The Second part)}

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.

\bd{athletes}\ \
{\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.
\ed

\bd{time-int.dfn}\ \
Set ($\epsilon$) and ($c$) to some desired values depending on the time
period of the set of classes, now, the above {\it effective}
$t$-neurons can be defined as follows:

\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
\ed

\bd{consts-dfn}\ \
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 in definition (\ref{dfn-5})
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.

\ed

\subsection{Constraints}

The study case dealt with here is {\it course scheduling} at the university level
for both graduate and undergraduate students. Given the various
student and course section information, a host of constraints need to
be satisfied for any successful scheduling of this nature. In addition
to those previously stated in this report, here is a sample ({\it hard and soft}):
\begin{enumerate}
\item Section must be open. ({\it hard})
\item Student must meet eligibility criteria for the section. ({\it
hard})
\item Section meeting time may not conflict with other ``enrolled''
sections. ({\it hard}).
\item An enrollment balance by gender should be obtained in Honors
seminar sections. ({\it soft})
\item For each student, balance $MWF$ and $TuTh$ sections. ({\it hard})
\item Do not enroll athletes in courses that conflict with their sport
practice time. (may be {\it hard} or {\it soft}, depending on the sport)
\item Balance enrollment in multi-section courses. ({\it soft})
\item Enroll students in ALL parts of a course when associated sections 
(i.e. recitation, labs) are required. ({\it hard})
\item Optimize 1st choice course enrollments by attempting all 1st
choice courses for all students before attempting 2nd and 3rd
choices. ({\it soft}) 

%%%%%%% \item Do not enroll students in sections when 
%%%%%%% {\it PRINT} is set to {\it NO}. 

\item A student with mapping status of {\it Req} should occupy a class
seat over a student with a status of priority ($\chi$). ({\it hard})
\item Class size/capacity to be observed. ({\it hard})
\item Distance/building preference or department-to-building distance
({\it soft}).
\item Various preferences of students or/and professors have to be
considered. ({\it soft}).
\end{enumerate}

Other constraints that may need to be taken into account during the
scheduling process:
\begin{itemize}
\item It is not the case, in general, that all classrooms are
available for all subjects. In reality, many subjects require a
special purpose rooms.
\item Some classes are taught for two or more hours in a row.
\item Group formation: for some subjects the classes are broken up
into option groups temporarily forming new classes.
\item The schedule should have as few ``holes'' as possible; that is
classes should be ``glued'' together.
\item Student with a ``required'' mapping status should occupy a class
seat over a student with other status.
\item Two courses of the same topic are scheduled on two different
working days (may be).
\item The predecessor-successor relation between the subjects/topics
is satisfied.
\item Classes may be designated as clashing with other classes because
of common students. 
\item The timetable may consist of an arbitrary number of periods per
day.
\item Some rooms may be declared unavailable during any time period.
\item Class subject combinations may request a preference for certain
periods.
\item Individual professor and class preferences may be specified for
specific rooms even if they are taken from a group or rooms.
\item Some classes may be grouped with others.
\item Daily limits must be specified to restrict professors from
taking more than a specified number of classes per day.
\item In order to schedule a number of classes at the same time, it is
possible to link class tuples together. In this way, the cost of a
move is the sum of the individual costs of each class tuple, however,
they are all moved at the same time.
\end{itemize}

\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\}$.

\bd{balance-dfn.1}
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
\ed

\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}.

\bd{athletic}

\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
\ed

\subsubsection{Eligibility for Enrollment}
\bd{The Eligibility Criteria}
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
\ed

\subsubsection{Course Pre-requisite and Subsequent Assignment}

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

\bd{pre-sub.dfn}\ \
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.
\ed

\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
non-Euclidean type CO problems such as 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

\subsection{Some Analysis of the MFA Performance on the TTS}

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, 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 these constraints are hard, but the majority are
classified to be soft. This later point can be compensated for by more
stringent annealing parameters, but at the cost of longer simulation
times. No attempts have been make to carry this out.
Also, we observe that the overall quality of the solution degrades as
the problem size increases. This degradation moves from soft
constraint violation to major violations of hard constraints, and
hence, illegal timetables.
From various correspondences with Andrew Gee, I'd like postulate
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 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 is no objective,
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 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
convergence to the optimal solutions for asymptotically slow
annealing, such results is not true for MFA algorithms, because MFA is
an approximation.

\section{Highlights of Mean Fields Annealing}

\subsection{Brief Background}

To get around the difficulties arising from Hopfield-Tank (HT)
algorithm, one way forward is to reduce the computational burden being
placed on the network:

{\em Constrain each neuron representing , for example, a
class-professor to be ``on'' only once (each class-professor to
take place only once during a certain time period.) by enforcing:} 

\be{normal}
{\displaystyle \sum_a} V_{ia} = 1\;\; \forall i, V_{ia} = 1/0\;;
\ee

{\em saying that event $i$ to take place during interval $a$ or $i$ is the
$a$-th event to occur in the schedule, rather than relying
on an energy penalty term to try to do this (as is done in HT)}.

\subsection{MFA variables: Ising to the Potts model}

\subsubsection{The Ising Case}

The Ising model describes a magnetic system in terms of binary spins
$s_i \in \{-1 , 1\}$. The system is governed by the energy function

\be{Ising}
E = - {J \over 2} {\displaystyle \sum_{i}}\, s_i s_{i+1}\;\;,
\ee

where nearest neighbors interact pairwise with a constant attractive
(if $J > 0$) coupling of strength $J$. The lowest energy state is
reached by iterative updating of 

\be{interact}
s_i = sgn[J(s_{i-1} + s_{i+1})]\;,
\ee

which leads to a state where all spins point in one of two possible
directions.

\subsubsection{The Potts Representation}

Restricting the state space of the neurons such that exactly {\it one
neuron} at every site is allowed to be on, can be written as: 

\be{norm.1}
{\displaystyle \sum_{a}} s_{ia} = 1\;\;\;\;\forall i\;.
\ee

Thus for every $i$, $s_{ia}$ should be 1 for only one value of $a$,
and $0$ for the rest. 

Choose a $k$-dimensional vector representation of the neuron states:
The $a^{th}$ possible state of the neuron vector $\vec{s}_i$ associated
with the node $i$ is represented by the unit vector in the direction
of $a$, so we end up with $K$ different states, 
$(1,0,\ldots,0)$, $(0,1,\ldots,0)$, $\ldots$, $(0,0,\ldots,1)$ are all
normalized and mutually orthogonal. The number of states available at
every node is thereby reduced from $2^K$ to $K$. Technically, we now have
a $K$-state Potts model.

The explicit enforcement of the following equality is referred to as
{\it neuron normalization}:

\be{norm.2}
{\displaystyle \sum_{a \in A}} V_{ia} = 1\;,
\ee

where $A$ is a set of states, in the MFA algorithm without any need
for a penalty term.

\subsection{Mapping Timetabling onto the Potts Neural Nets}

\subsubsection{The Relationship Between the Key Timetabling Entities}

Our problem consist of the following sets: teachers ($P$), classes
($C$), students ($S$), rooms ($R$), and time periods (${\cal I}$).
\vspace{0.1cm}
\begin{figure}[htb]
\begin{center}
\begin{picture}(230,230)
\psfig{figure=scheduling-net.eps,height=230pt,width=230pt}
\end{picture}
\end{center}
\caption{A one-way arrow indicates that the ``bottom'' entity depends
on the ``top'' one. A two-way arrow indicates a mutual dependency
between the two entities.}
\label{relations}
\end{figure}

\subsubsection{Neural Mapping}

Involve mapping the problem onto the Potts model.

\begin{itemize}
\item These multi-state neurons, forming the network, represent different
activities of the timetabling problem. 

\item Construct the constraints of the problem using these neurons.

\item Derive those MF equations controlling the dynamics of the net:

\begin{itemize}
\item Introduce the MF variables , for example, 

\be{rel.1}
V_{pq;xt} = <S_{pq;xt}>_T\;\;.
\ee

\item Start off with the partition function:
\be{part}
Z = {\displaystyle \sum_{[\vec{s}]}} e^{-E[\vec{s}] / T} \;.
\ee

\item Re-write $Z$ in terms of a multidimensional integral
over continuous vectors ($\vec{u}_i, \vec{v}_i$).

\item The MFT approximation to $<\vec{s}_i>$ is given by the value
$\vec{u}_i$ at a saddle-point equations:\\

\be{mft.1}
{{\partial E} \over {\partial v_{ia}}} = 0
\ee
and
\be{mft.2}
{{\partial E} \over {\partial u_{ia}}} = 0\;,
\ee
generating a set of self-consistency equations in the MF variables:

\be{mfv.1}
u_{ia} = 
- {1 \over T} {{\partial E[\vec{v}]} \over {\partial v_{ia}}}\;,
\ee
and
\be{mfv.2}
v_{ia} = {e^{u_{ia}} \over {{\displaystyle \sum_b} e^{u_{ib}}}}\;.
\ee

Where $a$ denotes the vector components.

For events (pq, xt), we get the following MF variables:

\be{mfv.3}
u_{pq,xt} = 
- {1 \over T} {{\partial E} \over {\partial v_{pq,xt}}} \;,
\ee
and
\be{mfv.4}
v_{pq,xt} = 
        {e^{u_{pq,xt}} \over {{\displaystyle \sum_{x' t'}}
e^{u_{pq,x't'}}}} \;.
\ee

From the later follows, that ${\displaystyle \sum_{a}} v_{ia} = 1\;\;$.
One can think of the mean field $v_{ia}$ as the probability for the
Potts neuron $i$ to be in state $a$.
\end{itemize}

\item Iterate $v_{pq, xt}$ equation with annealing, i.e. starting at a
high temperature, and successively lowering it in the course of the
process.  A phase transition is passed at $T = T_c$.

\item Estimate the critical temperature $T_c$, which sets the scale of
$T$, by expanding the above $v_{pq,xt}$ equation around the trivial
fix-point : at $T = 0$
\be{fixed}
v^{(0)}_{pq,xt} = {1 \over N_x N_t}\;.
\ee
\end{itemize}

\section{Timetabling and Simulated Annealing (SA)}

Simulated annealing has been widely used for tackling different
combinatorial optimization problems, particularly academic
scheduling. The results obtained depend heavily on the cooling
schedule used. We initially used the most commonly known and used
schedule, which is the geometric 
cooling schedule, but later tried the geometric re-heating
schedule based on cost, due to Abramson and colleagues \cite{2}. 

A comprehensive discussion of the theoretical and practical details of
simulated annealing is given in \cite{17}. It
suffices here to say that the elementary operation in the Metropolis
method for a combinatorial problem such as scheduling is
the generation of some new candidate configuration, which is then
automatically accepted if it lowers the solution cost ($C$), or accepted
with probability $exp(-\Delta C /(k T))$, where $k$ is the Boltzmann
constant and $T$ is the temperature, if it would increase the cost by
$\Delta C$. Thus the technique is essentially a generalization of the
local optimization strategy, where, at non-zero temperatures, thermal
excitations can facilitate escape from local minima.

\subsection{Some Advantages of SA}

\begin{enumerate}
\item Relatively easy to implement.
\item Applicable to almost any combinatorial optimization problem.
\item It usually provides reasonable solutions, depending on the
cooling schedule(s) used. 
\end{enumerate}

\subsection{The Algorithm}

The following is a sketch of the SA algorithm with a ``generic''
cooling schedule:
\vspace{0.5cm}
Let ${\cal N}_{s}$ be the neighborhood of configuration $s$. $k$ is
the Boltzmann constant $=$ 1.
\begin{enumerate}
\item Generate an initial solution (configuration) $s$ .
\item Compute cost of $s$ :  $C(s)$.
\item Compute initial temperature $T_0$.
\item Set the temperature $T\; =\; T_0$.
\item While \underline{stop criterion} is not satisfied do:
\begin{enumerate}
\item Repeat \underline{Markov chain length} (M) times:
\begin{enumerate}
\item Select a configuration $s^{'} \subset {\cal N}_{s}\;$.
\item Let $\Delta = C(s^{'}) - C(s)\;$.
\item If ($\Delta \leq 0$ $\{$downhill move$\}$), set $s = s^{'}\;$.
\item If ($\Delta > 0$ $\{$uphill move$\}$), set $s = s^{'}$ with a
probability $= e^{-\Delta /(k T)}\;$.
\end{enumerate}
\item Reduce (or update) temperature $T$.
\end{enumerate}
\item Return the solution  $s$.
\end{enumerate}
\vspace{0.25cm}
\hrule
\hrule
\vspace{0.5cm}

\subsection{Mapping of Timetabling Into Annealing Algorithm}

The most obvious mapping of timetabling (TTS) into the SA algorithm
involves the following constructs:

Let max(P) = maximum number of classes assigned to any professor.
\begin{enumerate}
\item a {\bf state} is a timetable containing the following sets:
\begin{itemize}
\item P: a set of professors.
\item C: a set of classes.
\item S: a set of students.
\item R: a set of classrooms.
\item ${\cal I}$: a set of time periods.
\end{itemize}
\item a {\bf cost} or ``energy'' E(P, C, S, R, ${\cal I}$) such that:
\begin{itemize}

\item E(P): is the cost of assigning more than max(P) classes to the
same professor, plus scheduling one or more courses that causes a conflict
in the professor's schedule.

\item E(C): is the cost of scheduling certain classes at/within the
same time period in violation of the exclusion constraint, for
example.

\item E(S): is the cost of having two or more classes conflict in
time; plus cost of having in the schedule one or more classes that
really don't meet the student's major, class requested, or course
requirements; plus the cost of not having the courses evenly spread
over the week, etc. 

\item E(R): is the cost resulting from assigning room(s) of the wrong
size and/or type to a certain class.

\item E(${\cal I}$): is the cost of having more or less time periods than
required, plus cost of an imbalanced class assignments (a certain
period will have more classes assigned to than others, etc.).
\end{itemize}

\item a {\bf swap} (or a move) is the exchange of one or more of the
following:

\begin{itemize}
\item Class swapping: class $c_i$ with class $c_j$ in the set $C$ with
respect to time periods ${\cal I}_i$ and ${\cal I}_j$, and/or with respect to
classroom $R_i$ and $R_j$, respectively.
\end{itemize}
\end{enumerate}

\subsection{The Annealing Schedules}

Two annealing schedules have been used for our experiments: geometric
cooling and an adaptive re-heating as a function of cost (RFC, for
short). This later schedule is due to Abramson\cite{2}. 

\subsubsection{Geometric Annealing Schedule}

The new temperature in the above schedule is usually calculated
according to:

\be{cooling}
T^{'}\; =\; \alpha\; T\;,
\ee

where $\alpha$ is typically in the range $0.95 \ldots 0.99$. 
This schedule is used as a base line for comparison of other
scheduling schemes.

{\bf \underline{Advantages}}:
\begin{itemize}
\item Have a good theoretical foundation.
\item Proven history of being good in practice.
\item Well understood by all researchers working on annealing.
\end{itemize}

{\bf \underline{Disadvantages}}:

The state of the system is not taken into account in any way, and
thus can not adapt the intensity of the search depending on the
complexity of the problem.

\subsubsection{Reheating as a Function of Cost (RFC) Schedule}

{\bf \underline{Some background}}:
\begin{itemize}
\item In combinatorial optimization problems, {\bf phase transitions} can be
observed as sub-parts of the problem are resolved. Lister
\cite{7} , in tackling the TSP, and Abramson
\cite{2} observed that at high temperatures the gross
structure of the solution is resolved, and at low temperatures the
minute details of the solution are resolved.

\item Applying a reheating type procedure, depending on the phase,
would allow the algorithm to spend more time in the low temperature
phases, thus reducing the total amount of time required to solve a
given problem.

\item To calculate the temperature at which a {\bf phase transition}
occurs, it is necessary to compute what is called the {\bf specific
heat} of the substance. A {\bf phase transition} occurs at a
temperature $T_{msh}$ when the specific heat is maximal, and this
triggers the change in the state ordering. 

\item If the best solution found to date is poor (i.e. has a high
energy or cost) then the super-structure may require
re-arrangement. This can be done by raising the temperature to a level
which is higher than the phase transition temperature
($T_{msh}$). Generally, the higher the best cost, the higher the
temperature which is required to escape the local minimum.

\item $T_{msh}$ is the temperature at which the maximum specific heat
occurs, or at which the system changes state.
\end{itemize}


%%Before summarizing this particular schedule, we would like to state a
%%very striking resemblance between it and MFA algorithm. Both rely on
%%a specific heat or temperature that is computed before hand and at
%%which a first-order phase transition occurs (see e.g.
%%\cite{Mouritsen-84} for a discussion of phase transitions). For MFA,
%%this specific heat is the $T_c$ and for the RFC schedule is denoted as
%%$T_{msp}$ at which the maximum specific heat occurs. This $T_{msp}$ is
%%computed as follows (see \cite{Abramson-94} for more details):

{\bf \underline{The schedule}}:

\begin{itemize}
\item Let $C_i$ to denote the cost of configuration $i$, and $k$ the
Boltzmann constant.

\item At temperature $T$, the acceptance probability is: 
\be{T.1}
Pr_i(T) = {e^{-C_i \over {kT}} \over \sum_j{e^{-C_j \over {kT}}}}\,\,.
\ee

\item For a set of configurations ${\cal N}$, the expected cost at
temperature $T$: 
\be{T.2}
<C(T)> = \sum_{i \in {\cal N}} C_i Pr_i(T)\,\,.
\ee

\item The expected square cost at $T$: 
\be{T.3}
<C^{2}(T)> = \sum_{i \in {\cal N}} C_i^2 Pr_i(T)\,\,.
\ee

\item The variance of the cost (or energy) at $T$ is:
\be{T.4}
\sigma^{2}(T) = <C^{2}(T)> - <C(T)>^{2}\,\,.
\ee

\item Using the above derivations, the specific heat: 
\be{T.5}
T_{sh}(i) = {\sigma^{2}(T) \over T^{2}}\,\,, \text{where}\,i \in
\{1,2,\ldots,n\,\}\,.
\ee

\item Now the maximum specific heat:
\be{T.5.1}
T_{msh} = max(T_{sh}(i))\,\,, \text{where}\,i = 1,2,\ldots,n\,.
\ee

\item The new temperature is set to:
\be{T.6}
T_{new} = P \cdot C_b(s) + T_{msh}\,\,.
\ee

$T_{new}$ is the reset temperature, $P$ is a constant, and
$C_b(s)$ is the current best cost of the given configuration $s$.
\end{itemize}

It follows from the Boltzmann distribution (see White \cite{41}
and lemma 4.1.10 of \cite{13}) that 
\be{T.7}
{d \over {dT}}E_T = {d \over {dT}}C(T) = {\sigma^{2}(T) \over T^{2}}
\ee

%%\vspace{0.25cm}
%%\begin{figure}[htb]
%%\centerline{\psfig{figure=Temperature.vs.Tmsp.ps,height=320pt,width=320pt,angle=-90}}
%%\caption{Temperature vs Specific Heat}
%%\label{graph.2}
%%\end{figure}
%%
%%\vspace{0.25cm}
%%\begin{figure}[htb]
%%\centerline{\psfig{figure=Temperature.vs.Acceptance.ps,height=320pt,width=320pt,angle=-90}}
%%\caption{Temperature vs Percentage of Acceptance}
%%\label{graph.3}
%%\end{figure}

\section{More on Simulated Annealing}
Simulated annealing is a robust technique and has been applied with
success to a wide variety of problems. However, to obtain good results
the annealing runs often require a great deal of computer time.
\cite{19} criticized much of the earlier research for only
comparing simulated annealing with the simpler local optimization
heuristics (i.e. those using r-opt style moves) and not with the
more competitive Lin-Kernighan (LK) heuristics. They undertook a
systematic comparison of the above three methods' performance on the
graph bisection problem, for low degree random and random geometric
graphs. 
When assessed using equivalent amounts of computer time, simulated
annealing emerged a clear winner over simple local optimization for
both types of graph. Against the LK-style algorithm however, the new
technique (SA) was only slightly better for random graphs and worse
for the random geometrical class. A similar investigation is currently
underway for the TSP \cite{21, 22}; preliminary results
\cite{20} suggest that, matched for run-time, simulated
annealing beats 3-opt (see Table.\ref{table.6}) 
but is worse than LK. 

In summary, simulated annealing appears competitive with
many of the best local optimization methods. Its major drawback
is the requirement for long run times. Furthermore, the inherently
serial nature of the approach in general precludes any obvious
parallel implementation, though several attempts have been made to
overcome this problem such as the recent work by Frost \cite{16}.

\section{Simulation Issues and Results}
Our study case is quite large in comparison to the scheduling case used
by Peterson and colleagues \cite{28}. The largest case they
have dealt with involve around 1000 students, 20 different possible
majors, and an overall periodic schedules (over weeks). In the case of
Abramson and colleagues \cite{2}, the problems they dealt
with are quite small and they state that problems involving more than
300 tuples\footnote{A tuple is a combination of a class, teacher, and
a classroom. The time periods are assumed to be given.} are quite
difficult to solve.

Our study case involve real scheduling data covering three semesters
here at Syracuse University. 
Those of first and second semester are roughly of the same size:
around 13653 students, 2500 of professors, staff, and teaching
assistants, 581 rooms (classrooms, seminar rooms and auditoriums), and 
3624 classes. With the exception of number of rooms, the third
semester data is much smaller. Along with this, we also have utilized 
the distance (in a form of a distance matrix) between various
buildings on main campus. Given the various building preferences of
the departments, we have constructed a ``finite state automaton'' in a
matrix format $[I,J]$ between all the academic schools and departments
$[I]$ and all the buildings involved in scheduling $[J]$. Using this
matrix in conjunction with the real -appropriately scaled- distance, a
final state matrix is derived and directly used in the scheduling
process. 

\subsection{On Results obtained via SA and MFA}

The main objective here is to generate a solution (i.e. a feasible
timetable) subject to all binding (or hard) constraints and to
minimize the cost incurred form the medium and soft constraints.


\begin{table}[htb]
\vspace{0.25cm}
\hrule
\hrule
\vspace{0.25cm}
The following table illustrates some mean error results of the MFA and
Simulated annealing using the RFC schedule, with and without
the rule-based system as a pre-processor. Each of these figures is for
100 runs of the same size (class, classrooms, students, etc.)
timetabling problem.
\vspace{0.25cm}
$$
\begin{tabular}{|l||c|c||c|c||}
\hline
                 &\multicolumn{2}{|c|}{Mean Field Annealing}&\multicolumn{2}{|c|}{Simulated Annealing} \\ \hline\hline
                 & Pre-Processor & No Pre-Processor & Pre-Processor   & No Pre-Processor \\ \hline\hline
Mean        &  2.6   &  14.24 &  0.32  & 7.1   \\
Percentage  &        &        &        &       \\       
Error       &        &        &        &       \\ \hline
\end{tabular}
$$
\caption{Performance figures for MFA and SA in the presence and
absence of pre-processor (a rule-based knowledge system).}
\label{table.1}
\vspace{0.25cm}
\hrule
\hrule
\end{table}

\begin{table}[htb]
\vspace{0.25cm}
$$
\begin{tabular}{|l||c|c|c|c||}
\hline
            & Total Courses & Scheduled & Not Scheduled & \% Scheduled\\ \hline\hline
First Term  &     4400      &   2992    &      1408     &     68.00   \\ \hline
Second Term &     4082      &   2881    &      1201     &     70.57   \\ \hline
Third Term  &     709       &   569     &      140      &     80.25   \\ \hline
\end{tabular}
$$
\caption{Performance figures for three terms using a {\it rule-based
system} only. Both class and students constraints are used.}
\label{table.2}
\vspace{0.25cm}
\hrule
\hrule
\end{table}

\begin{table}[htb]
\vspace{0.25cm}
$$
\begin{tabular}{|l||c|c|c|c||}
\hline
            & Total Courses & Scheduled & Not Scheduled & \% Scheduled\\ \hline\hline
First Term  &     4400      &   3168    &      1232     &     72.00   \\ \hline
Second Term &     4082      &   3189    &      893      &     78.12   \\ \hline
Third Term  &     709       &   631     &      78       &     89.00   \\ \hline
\end{tabular}
$$
\caption{Performance figures using a {\em rule-based system (or expert
system)}. The only type of constraints used are those for classes.} 
\label{table.3}
\vspace{0.25cm}
\hrule
\hrule
\end{table}

\begin{table}[htb]
\vspace{0.25cm}
$$
\begin{tabular}{|l||c|c|c||}
\hline
                & Scheduled & Not Scheduled & \% Scheduled\\ \hline\hline
SA (geometric)  &   2627    &      1773     &     59.70   \\ \hline
SA (cost-based) &   2914    &      1486     &     66.22   \\ \hline
Mean Field      &   2151    &      2249     &     48.88   \\ \hline
\end{tabular}
$$
\caption{Performance figures for three algorithms with no
pre-processor used. Both students and class constraints are used and
the data is that of the first term. \underline{NOTE}: Annealing is
carried out on space AND time. Also, all shown figures are averaged
over 10 runs.}
\label{table.4}
\vspace{0.25cm}
\hrule
\hrule
\end{table}

\begin{table}[htb]
\vspace{0.25cm}
$$
\begin{tabular}{|l||c|c|c||}
\hline
                & Scheduled & Not Scheduled & \% Scheduled\\ \hline\hline
SA (geometric)  &   4280    &      120      &     97.27   \\ \hline
SA (cost-based) &   4333    &      67       &     98.47   \\ \hline
Mean Field      &   3985    &      415      &     90.56   \\ \hline
\end{tabular}
$$
\caption{Performance figures for three algorithms using an expert
system as a pre-processor. Both students and class constraints are
used and the data is that of the first term. \underline{NOTE}:
Annealing is carried out on space AND time. Also, all shown figures
are averaged over 10 runs.}
\label{table.5}
\vspace{0.25cm}
\hrule
\hrule
\end{table}

\begin{table}[ht]
The following table illustrates some mean error results of the four
models on  1000 runs of a Non-Euclidean 10-city TSP (communicated by
Andrew Gee). The figures shown here for neural nets are quite
discouraging.
\vspace{0.25cm}
$$
\begin{tabular}{|l||c||c||c||c||}
\hline
            &{Hopfield Net}&{MFA}&{Simulated Annealing}&{3-opt search}\\ \hline\hline
Mean        &  13.5  &  9.73  &  0.29  & 0.84  \\
Percentage  &        &        &        &       \\       
Error       &        &        &        &       \\ \hline
\end{tabular}
$$
\caption{Performance figures for Hopfield-net, MFA, SA, and 3-opt on a
non-Euclidean TSP case.}
\label{table.6}
\vspace{0.25cm}
\hrule
\hrule
\end{table}

%%--------------------------------------------------------------------
%% Need to include plots: temperature.vs.cost,
%% professor_cost.vs.timesteps, temperature.vs.delta_cost, 
%% roomsizecost.vs.timesteps, class.1.vs.rmsizecost, class.2.vs.rmsizecost,
%% diff.in.cost_class1-class2(pair of swapped classes),
%% class1_enroll.vs.rmsizecost,  class2_enroll.vs.rmsizecost,
%% class1_capacity.vs.rmsizecost, class2_capacity.vs.rmsizecost
%%

\subsection{On Results obtained via MFA}

To summarize, the results obtained from mean field annealing range
from not quite satisfactory to totally unacceptable. Not quite
satisfactory in the sense that when dealing with either hard or soft
constraints separately (but not at the same time), the resulting
timetables do satisfy the given constraints, but depends on how many
of these constraints we're dealing with. As the number of constraints
increases, the violations start appearing till it is no longer viable
to continue with the procedure. 

When running the method on a combined hard and soft constraints,
violations start occurring at a higher frequency and many involves the
hard constraints, clearly, yielding illegal timetables. 


\section{Summary and Conclusion}
This part outlines the basic method and the scheduling constraints.
An academic scheduling problem that is based upon equality constraints
is considered. 
Two methods are used to apply our data and produce the target
schedules. This target need to be a valid timetable of classes,
professors, and students with no violations at all of the hard
constraints and preferrabaly satisfying most , if not all, of the soft
constraints. Using simulated annealing with an adaptive re-heating
schedule based on cost and with no pre-processing, over 60 out of each
100 timetables produced, do satisfy all of the hard constraints and
most of the soft constraints used. On the other hand, mean field
annealing produces less than 50 legal timetables out of each 100
timetables produced. The rule-based system which is used as
pre-processor for each method contributed greatly to the improvement
of the simulated annealing results as well as those of mean field
annealing. The use of pre-processing with MFA improved the results to
around 90 percent; that is 90 out of 100 schedules do satisfy all of
the hard constraints. While the use of pre-processing with SA yielded
a 98 percent accuracy; that is 98 out of 100 schedules are legal
(satisfying all the hard constraints and most of medium and soft
ones).

In lights of these results, the figures reported in
Table-\ref{table.1} seem to be quite consistent with the figures of
Table-\ref{table.6} of the two methods, simulated annealing and MFA on
the class (represented by the TSP) of non-Euclidean CO problems, which
the scheduling problem belongs to. 

We plan to carry out more experiments using the techniques reported in
this paper, perhaps using an improved pre-processor system.

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

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

\bibitem{2}
Abramson, Dang, and Krishnamoorthy,
{\it An Empirical Study of Simulated Annealing Cooling Schedules},
Griffith Univ. report, Nathan, Qld, Aus. 1994.

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

\bibitem{4}
J. Thompson, 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{5}
K. Gosselin and M. Truchon, ``Allocation of Classrooms by Linear
Programming'', Journal of Operational Research Society 37, 561 (1986).

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

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

\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{10}
Ladkin, P.B., and Maddux, R.D.,
``On binary constraint problems'', {\it Journal of the ACM} 41 (1994),
pp. 435-469.

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

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

\bibitem{13}
Sorkin, Gregory,
{\it Theory and Practice of Simulated Annealing on Special Energy
Landscapes}. PhD Thesis, UC-Berkeley, Dept of EE and CS, 1992.

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

\bibitem{15}
Shrijver, A,
{\it A Course in Combinatorial Optimization}, to be published in 1995.

\bibitem{16}
Forest, Richard,
{\it Ensemble of Simulated Annealing}, SDSC, California, 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{19}
Johnson, Aragon, McGeoch, and Schevon,
Optimization by Simulated Annealing: an Experimental Evaluation, Part
I (Graph Partitioning). OR 37, 865-892 (1989).

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

\bibitem{21}
Johnson, Aragon, McGeoch, and Schevon,
Optimization by Simulated Annealing: an Experimental Evaluation, Part
III (TSP). Work in progress.

\bibitem{22}
David Johnson,
Data Structures for the TSP, to appear in the Jr. of Algorithms.

\bibitem{23}
F. Praparata, M. Shamos, 
{\it Computational Geometry}, Springer-Verlag, 1988.

\bibitem{24}
Minsky, Marvin and S Papert, {\it Perceptron} , 1988 edition.

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

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

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

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

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

\bibitem{30}
Peterson, 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 , second half of 1993.

\bibitem{31}
Paul Stolorz, 
``Abusing Statistical Mechanics to do Adaptive Learning and Combinatorial
Optimization'', Preprint, Santa Fe Institute, Dec 1991.

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

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

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

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

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

\bibitem{37}
A. Colorni, 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{38}
Alain Hertz, ``Finding a feasible course schedule using Tabu search'',
Discrete Applied Mathematics 35, 255-270, 1992.

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

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

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

\bibitem{42}
Martin Simmen,
PhD Thesis, Univ. of Edinburgh, 1993.

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

\bibitem{44}
Martin Simmen, private communication.

\bibitem{45}
Andrew Gee, private communication.

\bibitem{46}
Carsten Peterson, private communications.

\bibitem{47}
Mike Carter, private communications.

\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}

