課程設(shè)計交巡警服務(wù)平臺的設(shè)置與調(diào)度

1、交巡警服務(wù)平臺的設(shè)置與調(diào)度摘 要對提供的數(shù)據(jù)進(jìn)行了預(yù)處理后，利用djist算法求出任意兩個路口節(jié)點之間的最短路程，程序見附錄2，結(jié)果見“a區(qū)節(jié)點間最小路的值.xls”。交巡警服務(wù)平臺分配管轄范圍的原則是（1）所管轄的節(jié)點到該交巡警服務(wù)平臺的最近路程應(yīng)盡量小于3000m,（2）各交巡警服務(wù)平臺管轄范圍內(nèi)的所有節(jié)點發(fā)案率之和應(yīng)盡量接近平均發(fā)案率6.225。由此建立多目標(biāo)規(guī)劃模型。對該模型利用matlab編程（見附錄1）求出滿意解如圖1所示。對進(jìn)出該區(qū)的13條交通要道快速封鎖問題，我們把參與封堵各交巡警服務(wù)平臺到各出口要道最短路程之和的最小化作為第一目標(biāo)，把參與封堵各交巡警服務(wù)平臺到各出口要道的最長

2、路程最小化作為第二目標(biāo)，建立多目標(biāo)規(guī)劃模型，利用matlab編程（見附錄4）求解，結(jié)果如表1和圖2所示。根據(jù)a區(qū)各交巡警服務(wù)平臺的出警時間方差和工作量方差，考慮在服務(wù)平臺15附近的a區(qū)出口節(jié)點29增設(shè)一個服務(wù)平臺，在服務(wù)平臺7附近的路口節(jié)點48增設(shè)一個服務(wù)平臺，在路口節(jié)點91增設(shè)一個服務(wù)平臺。增設(shè)3個服務(wù)平臺后的a區(qū)巡警服務(wù)平臺分配管轄范圍如圖3所示，合理性有較大改善。對全市六區(qū)交巡警服務(wù)平臺設(shè)置合理性評價如a區(qū)，主要要避免服務(wù)平臺的出警時間過長，不能及時趕到突發(fā)事件地點，并且使各服務(wù)平臺的工作量盡量均衡。我們計算了b區(qū)的出警時間方差和工作量方差如表3所示（程序見附錄6）。建議在99、94、9

3、8處進(jìn)行增設(shè)服務(wù)平臺，另外也可以考慮在95處增加警力。其他各區(qū)的處理方法同b區(qū)。該市p點發(fā)生了重大刑事案件，快速圍堵嫌疑犯的調(diào)度問題，是一個動態(tài)規(guī)劃。我們給出了嫌疑犯逃跑各階段的路口集合，可能的逃跑路線如圖4所示。用matlab編程（程序見附錄7），得到最佳圍堵方案如表4。關(guān)鍵字：動態(tài)規(guī)劃 多目標(biāo)規(guī)劃 最短路徑算法 圖論 0-1規(guī)劃一、 問題的重述“有困難找警察”，是家喻戶曉的一句流行語。警察肩負(fù)著刑事執(zhí)法、治安管理、交通管理、服務(wù)群眾四大職能。為了更有效地貫徹實施這些職能，需要在市區(qū)的一些交通要道和重要部位設(shè)置交巡警服務(wù)平臺。每個交巡警服務(wù)平臺的職能和警力配備基本相同。由于警務(wù)資源是有限的，

4、如何根據(jù)城市的實際情況與需求合理地設(shè)置交巡警服務(wù)平臺、分配各平臺的管轄范圍、調(diào)度警務(wù)資源是警務(wù)部門面臨的一個實際課題。試就某市設(shè)置交巡警服務(wù)平臺的相關(guān)情況，建立數(shù)學(xué)模型分析研究下面的問題：（1）根據(jù)給出的相關(guān)數(shù)據(jù)信息，為各交巡警服務(wù)平臺分配管轄范圍，使其在所管轄的范圍內(nèi)出現(xiàn)突發(fā)事件時，盡量能在3分鐘內(nèi)有交巡警（警車的時速為60km/h）到達(dá)事發(fā)地。對于重大突發(fā)事件，需要調(diào)度全區(qū)20個交巡警服務(wù)平臺的警力資源，對進(jìn)出該區(qū)的13條交通要道實現(xiàn)快速全封鎖。實際中一個平臺的警力最多封鎖一個路口，給出該區(qū)交巡警服務(wù)平臺警力合理的調(diào)度方案。根據(jù)現(xiàn)有交巡警服務(wù)平臺的工作量不均衡和有些地方出警時間過長的實際情

5、況，擬在該區(qū)內(nèi)再增加2至5個平臺，請確定需要增加平臺的具體個數(shù)和位置。（2）針對全市（主城六區(qū)a，b，c，d，e，f）的具體情況，按照設(shè)置交巡警服務(wù)平臺的原則和任務(wù)，分析研究該市現(xiàn)有交巡警服務(wù)平臺設(shè)置方案的合理性。如果有明顯不合理，請給出解決方案。如果該市地點某點處發(fā)生了重大刑事案件，在案發(fā)3分鐘后接到報警，犯罪嫌疑人已駕車逃跑。為了快速搜捕嫌疑犯，給出調(diào)度全市交巡警服務(wù)平臺警力資源的最佳圍堵方案。二、 問題的分析1、交巡警服務(wù)平臺分配管轄范圍的基本原則應(yīng)該是（1）所管轄的范圍內(nèi)出現(xiàn)突發(fā)事件時，應(yīng)盡量能在3分鐘內(nèi)有交巡警到達(dá)事發(fā)地，其實就是所管轄的節(jié)點到該交巡警服務(wù)平臺的最近路程應(yīng)盡量小于30

6、00m,（2）各交巡警服務(wù)平臺的工作量應(yīng)盡量均衡，就是各交巡警服務(wù)平臺管轄范圍內(nèi)的所有節(jié)點發(fā)案率之和應(yīng)盡量接近平均的發(fā)案率。由此可見該問題是一個多目標(biāo)規(guī)劃問題。為了求解問題，我們首先對提供的數(shù)據(jù)進(jìn)行一些預(yù)處理，再利用djist算法求出任意兩個路口節(jié)點之間的最短路程，程序見附錄二，結(jié)果見“a區(qū)節(jié)點間最小路的值.xls”。由于題目提供的節(jié)點坐標(biāo)單位為mm，因此我們將警車3分鐘走過的路程根據(jù)比例化為30mm。如該區(qū)發(fā)生重大事件，交巡警對進(jìn)出該區(qū)的13條交通要道要實現(xiàn)快速全封鎖，即是求各交巡警服務(wù)平臺到13條交通要道最小路徑的總體最小化，以滿足盡量實現(xiàn)快速封鎖的目的。也是一個多目標(biāo)規(guī)劃問題。增設(shè)交巡警

7、服務(wù)平臺的原則是：（1）對有些出警時間過長的地方，可以考慮增設(shè)交巡警服務(wù)平臺，以期盡量能在3分鐘內(nèi)有交巡警到達(dá)事發(fā)地。（2）有些工作量過大的交巡警服務(wù)平臺附近，可以考慮增設(shè)交巡警服務(wù)平臺，以期工作量應(yīng)盡量均衡。2、對該市現(xiàn)有交巡警服務(wù)平臺設(shè)置方案合理性的分析，主要是看目前的方案中是否存在有些服務(wù)平臺工作量過大和出警時間過長的情況，如果有明顯不合理，應(yīng)考慮增設(shè)服務(wù)平臺或?qū)ぷ髁窟^大的服務(wù)平臺增加警力。該市p點發(fā)生了重大刑事案件，在案發(fā)3分鐘后接到報警，犯罪嫌疑人已駕車逃跑。為了快速圍堵嫌疑犯，給出最佳調(diào)度方案，是一個動態(tài)規(guī)劃問題。 三、 模型假設(shè)1、任意相鄰兩個路口之間的道路都是直線。2、增設(shè)服

8、務(wù)平臺的時，其他服務(wù)平臺地點不變。3、不考慮交通堵塞，犯罪嫌疑人的車速與交巡警的車速相同，均為常量。4、各交巡警服務(wù)平臺的警力配備基本一致。四、 符號說明：表示從第i個服務(wù)平臺到第j個路口節(jié)點的最短路程。：表示從第i個服務(wù)平臺到第j個路口節(jié)點最短路程的正偏差。：表示第j個路口節(jié)點是否屬于第i個服務(wù)平臺，是為1否為0。：表示第i個服務(wù)平臺是否負(fù)責(zé)圍堵第j個交通要道，是為1否為0。：表示第j個路口節(jié)點的發(fā)案率。：a區(qū)平均每個服務(wù)平臺的發(fā)案率 ：表示從第k個服務(wù)平臺到第j個交通要道的最短路程。：表示第k個服務(wù)平臺是否負(fù)責(zé)封堵第j個交通要道，是為1否為0。五、 模型的建立和求解5.1 a區(qū)交巡警服務(wù)平

9、臺分配管轄范圍的優(yōu)化模型根據(jù)問題分析，交巡警服務(wù)平臺分配管轄范圍的原則主要是：（1）所管轄的節(jié)點到該交巡警服務(wù)平臺的最近路程應(yīng)盡量小于30mm，所以將總的正偏差最小化作為第一目標(biāo)。即最小化 其中：=（2）各交巡警服務(wù)平臺的工作量應(yīng)盡量均衡，因此將服務(wù)平臺所管轄的路口節(jié)點的總案發(fā)率的方差最小化作為第二目標(biāo)。即最小化建立如下的多目標(biāo)規(guī)劃模型：對該模型利用matlab編程（見附錄1）求出滿意解：如圖1所示：圖1 a區(qū)交巡警服務(wù)平臺分配管轄范圍示意圖5.2 交巡警服務(wù)平臺警力合理的調(diào)度模型對于該區(qū)發(fā)生重大事件，交巡警對13個交通要道快速全封鎖，根據(jù)問題分析，我們把參與封堵各交巡警服務(wù)平臺到各出口要道最

10、短路程的總和最小化作為第一目標(biāo)，把參與封堵各交巡警服務(wù)平臺到各出口要道的最長路程最小化作為第二目標(biāo)，建立多目標(biāo)規(guī)劃，模型如下： 其中 利用matlab編程（見附錄4）求解，結(jié)果如表1：表1 a區(qū)交巡警服務(wù)平臺對13條交通要道全封鎖方案出口節(jié)點標(biāo)號12222362483029281421241638服務(wù)平臺標(biāo)號12111345871514910162 圖2為示意圖，箭頭起點為服務(wù)平臺，終點為封鎖的出口節(jié)點。圖2 a區(qū)交巡警服務(wù)平臺對13條交通要道全封鎖方案示意圖5.3 a區(qū)交巡警服務(wù)平臺增設(shè)方案根據(jù)5.1交巡警管轄區(qū)域的劃分，通過計算得出a區(qū)各交巡警服務(wù)平臺的出警時間方差和工作量方差（見附錄5）

11、，結(jié)果如下。表2 a區(qū)各交巡警服務(wù)平臺的出警時間方差和工作量方差出警時間方差服務(wù)平臺標(biāo)號工作量方差服務(wù)平臺標(biāo)號010203060409050100601108012090140100160110190120.0380130.03150140.1440171.6310184.7330194.731816.47164.732036.15205.181346.5426.635141.66710.73171036.181522.807從表1中發(fā)現(xiàn)服務(wù)平臺15出警時間過長，服務(wù)平臺7出警時間過長且任務(wù)過重，服務(wù)平臺20出警時間較長且任務(wù)較重，因此考慮在服務(wù)平臺15附近的a區(qū)出口節(jié)點29增設(shè)一個服務(wù)平臺，

12、在服務(wù)平臺7附近的路口節(jié)點48增設(shè)一個服務(wù)平臺，在18和20之間的路口節(jié)點91增設(shè)一個服務(wù)平臺。增設(shè)三個服務(wù)平臺后的a區(qū)巡警服務(wù)平臺分配管轄范圍如圖3所示。圖3 a區(qū)增設(shè)后的巡警服務(wù)平臺分配管轄范圍示意圖5.4 全市交巡警服務(wù)平臺設(shè)置合理性評價及改進(jìn)方案對全市六區(qū)交巡警服務(wù)平臺設(shè)置合理性評價的依據(jù)如5.1，主要要避免服務(wù)平臺的出警時間過長，不能及時趕到突發(fā)事件地點，并且使各服務(wù)平臺的工作量盡量均衡。我們計算了b區(qū)的出警時間方差和工作量方差（程序見附錄6）。表3 b區(qū)交巡警服務(wù)平臺設(shè)置評價指標(biāo)出警時間方差服務(wù)平臺標(biāo)號工作量方差服務(wù)平臺標(biāo)號09309309509709709909801000100

13、0.25968.15964.419421.48945.2995229.569914.4498從表3中可以看出服務(wù)平臺99的出警時間過長，服務(wù)平臺94和96的出警時間較長，服務(wù)平臺98的工作量很大，服務(wù)平臺95和94的工作量也較大。因此可以采用5.3的方法在99、94、98處進(jìn)行增設(shè)服務(wù)平臺，另外也可以考慮在95處增加警力。其他各區(qū)的處理方法同b區(qū)。5.5 全市交巡警服務(wù)平臺警力資源的最佳圍堵嫌疑犯方案對于該市p點（第32個節(jié)點）處發(fā)生了重大刑事案件，在案發(fā)3分鐘后接到報警，犯罪嫌疑人已駕車逃跑。為了快速搜捕嫌疑犯，調(diào)度全市交巡警服務(wù)平臺警力資源的最佳圍堵方案如下：案發(fā)點為，以為始點的道路之終點

14、為，所有點的集合為稱之為狀態(tài)，為第一步可能的圍堵的路口集合，同理（集合元素的排序按距p的最短路程升序排列）稱之為狀態(tài)為第步可能的圍堵路口集合狀態(tài)轉(zhuǎn)移規(guī)律：1、在第步可能圍堵路口集合中，按升序依次判斷犯罪嫌疑人若從此路口逃跑是否能夠及時堵住（假設(shè)只要交巡警先到路口設(shè)卡，一定能夠堵住，到達(dá)時立即可以實施圍堵，設(shè)卡時間不計）。如果能堵住則繼續(xù)判斷下一個點，如果堵不住，則以該點為始點，向外擴展一步，即以該點為始點在逃逸方向上找終點，將若干個終點代替該點，得到狀態(tài)的圍堵點集合。2、如果某狀態(tài)中所有路口都能堵住，則可輸出最佳調(diào)度方案。通過編程計算第一狀態(tài)的路口集合為7，31，33，第二狀態(tài)的路口集合為8，

15、 34，37，47，15，30，第三狀態(tài)的路口集合為5，6，9，10，28，29，36，46，48，237，第四狀態(tài)的路口集合為16，26，35，39，49，45，49，50, 55，59，61，235，236，238，239，247，370，371，嫌疑犯從案發(fā)地點p可能逃跑的路線如圖4所示，不同的線型對應(yīng)著不同狀態(tài)下可能逃跑的路線及到達(dá)的路口節(jié)點。按照上述思路利用matlab編程計算的結(jié)果如表4（程序見附錄7）：表4 全市交巡警服務(wù)平臺圍堵嫌疑犯方案圍堵節(jié)點標(biāo)號569102829364648237服務(wù)平臺標(biāo)號56910151711687173p圖4 嫌疑犯從案發(fā)地點p可能逃跑的路線圖六、

16、模型的評價與推廣總的來說，模型與算法的思路清晰，其結(jié)果可幫助交巡警的調(diào)度和最佳行進(jìn)路徑的選擇，減少了出警時間，可操作性強，有一定的實際應(yīng)用價值，模型與算法也有一定的局限性和不足，如個別地區(qū)的出警時間超過3分鐘，也沒有考慮道路的暢通性等不定因素對出警時間的影響等.本模型與算法較好的解決了交巡警平臺的最優(yōu)選址問題，當(dāng)事故發(fā)生時，交巡警可以第一時間到達(dá)事發(fā)地點，有效的改善了交巡警在執(zhí)行任務(wù)中的效率，并且該模型與算法也可運用到其他應(yīng)急問題中，比如關(guān)于消防救援問題、重大生產(chǎn)安全事故應(yīng)急救援問題等.七、 參考文獻(xiàn)【1】 趙靜,但琦等編著，數(shù)學(xué)建模與數(shù)學(xué)實驗（第3版）， 北京，高等教育出版社，2008?！?

17、】 云舟工作室編著，matlab6數(shù)學(xué)建模基礎(chǔ)教程，北京，人民郵電出版社，2001。【3】 吳翊，吳孟達(dá)，成禮智編著，數(shù)學(xué)建模的理論與實踐，長沙，國防科技大學(xué)出版社，1999?！?】 譚永基，俞文吡編著，數(shù)學(xué)模型，上海，復(fù)旦大學(xué)出版社，1997。【5】 王海英等編著 圖論算法及其matlab實現(xiàn) 北京航空航天大學(xué)出版社 2010【6】 王桂平，王衍，任嘉辰編著 圖論算法理論、實現(xiàn)及應(yīng)用 北京大學(xué)出版社 2011附件附錄1%a區(qū)交巡警服務(wù)平臺分配管轄范圍求解clear all;mln = xlsread(d:b,2,a2:b138); %matrix of line between nodesm

18、n = xlsread(d:b,1,a2:e93); %matrix of nodes% create 0-1 matrix of line beteen nodesmln01=zeros(92,92);mlnsize=size(mln);mlnr=mlnsize(1); %row number of mlnfor i=1:mlnr if mln(i,1)=92 & mln(i,2)=92 % only deal with the data in the part a mln01(mln(i,1),mln(i,2)=1; endend% solve the path between the n

19、odesmlnpath=zeros(92,92);m=10000; %set m with a big numberfor i=1:92 for j=1:92 if mln01(i,j)=1 mlnpath(i,j)=sqrt(mn(i,2)-mn(j,2)2+(mn(i,3)-mn(j,3)2); end % set infinite number to the nodes with not connection endendsp=shortpath(mlnpath);for h=1:1000 randservicenode=randperm(20); handlenode=zeros(92

20、,1); %storage the depart result,0 note not servicenode=1:20; for i=1:20 startsn=servicenode(randservicenode(i); handlenode(startsn)=startsn; totalwork=0; clear totalpath startnodea totalpath(1)=0; startnodea(1,:)=startsn,0; startnodecount=1; while totalwork=10 startnode=startnodea(startnodecount,1);

21、 for j=1:mlnr if mln(j,1)=startnode & handlenode(mln(j,2)=0 tempm=startnodea(:,1); tempc=find(tempm=startnode); tempp=mlnpath(startnode,mln(j,2)+totalpath(tempc); if tempp30 totalpath(length(totalpath)+1)=mlnpath(startnode,mln(j,2)+totalpath(tempc); startnodea(size(startnodea,1)+1,:)=mln(j,2),startn

22、ode; handlenode(mln(j,2)=startsn; totalwork=totalwork+mn(mln(j,2),5); end end end startnodecount=startnodecount+1; if size(startnodea,1)sp(uanode,j) minservice=j; mindistance=sp(uanode,j); end end handlenode(uanode)=minservice; uam(i,2)=minservice; uam(i,3)=mindistance; end totalwork=0; largem=0; fo

23、r i=1:20 fservices=find(handlenode=i); work=0; for j=1:length(fservices) work=work+mn(fservices(j),5); end if work6.225 totalwork=totalwork+(work-6.225)2; largem=largem+1; end end finaltotalwork(h)=totalwork; finalhandlenodeh=handlenode; finallargem(h)=largem;enda,v=min(finaltotalwork)finallargem(v)

24、附錄2%求最短路徑函數(shù)function sp=shortpath(a)m=10000;a=a+a;for i=1:length(a) for j=1:length(a) if a(i,j)=0 a(i,j)=m; end endendfor i=1:length(a) pb(1:length(a)=0; pb(i)=1; d(1:length(a)=m; d(i)=0; temp=i; while sum(pb)length(a) tb=find(pb=0); d(tb)=min(d(tb),d(temp)+a(temp,tb); tmpb=find(d(tb)=min(d(tb); temp

25、=tb(tmpb(1); pb(temp)=1; end sp(i,:)=d;end附錄3%尋找a區(qū)出口節(jié)點clear all;mln = xlsread(d:b,2,a2:b928); %matrix of line between nodesmn = xlsread(d:b,1,a2:e93); %matrix of nodes% create 0-1 matrix of line beteen nodesmln01=zeros(92,92);mlnsize=size(mln);mlnr=mlnsize(1); %row number of mlnfor i=1:mlnr if mln(i

26、,1)=92 & mln(i,2)=92 % only deal with the data in the part a mln01(mln(i,1),mln(i,2)=1; endend% solve the path between the nodesmlnpath=zeros(92,92);m=10000; %set m with a big numberfor i=1:92 for j=1:92 if mln01(i,j)=1 mlnpath(i,j)=sqrt(mn(i,2)-mn(j,2)2+(mn(i,3)-mn(j,3)2); end % set infinite number

27、 to the nodes with not connection endendsp=shortpath(mlnpath);j=1;crossnode=zeros(13,1);for i=1:mlnr if (mln(i,1)92) | (mln(i,1)92 & mln(i,2)=92) if mln(i,1)=92 if isempty(find(crossnode=mln(i,1) crossnode(j)=mln(i,1); j=j+1; end else if isempty(find(crossnode=mln(i,2) crossnode(j)=mln(i,2); j=j+1;

28、end end end end附錄4%尋找封鎖交通要道的程序clear all;mln = xlsread(d:b,2,a2:b143); %matrix of line between nodesmn = xlsread(d:b,1,a2:e93); %matrix of nodes% create 0-1 matrix of line beteen nodesmln01=zeros(92,92);mlnsize=size(mln);mlnr=mlnsize(1); %row number of mlnfor i=1:mlnr if mln(i,1)=92 & mln(i,2)=92 % o

29、nly deal with the data in the part a mln01(mln(i,1),mln(i,2)=1; endend% solve the path between the nodesmlnpath=zeros(92,92);m=10000; %set m with a big numberfor i=1:92 for j=1:92 if mln01(i,j)=1 mlnpath(i,j)=sqrt(mn(i,2)-mn(j,2)2+(mn(i,3)-mn(j,3)2); end % set infinite number to the nodes with not c

30、onnection endendsp=shortpath(mlnpath);load out13.matsp13=zeros(13,20);for i=1:length(crossnode) tempsp13(i,:)=sp(crossnode(i),1:20); endsave shortpath.mat sp sp13for i=1:13 for j=1:20 if tempsp13(i,j)=0 tempsp13(i,j)=m; end endendfor h=1:1000 randoutnode=randperm(13); sp13=tempsp13; clear serviceton

31、ode servicetonode=zeros(13,3); % the first column store the out node,second store service node,third store distance servicetonode(:,1)=crossnode; % for i=1:13% if crossnode(randoutnode(i)=20% servicetonode(randoutnode(i),2)=crossnode(randoutnode(i);% servicetonode(randoutnode(i),3)=0; % end% end for

32、 i=1:13 if servicetonode(randoutnode(i),2)=0 a,b=min(sp13(randoutnode(i),:),2); tempm=servicetonode(:,2); while isempty(find(tempm=b) sp13(randoutnode(i),b)=m; a,b=min(sp13(randoutnode(i),:),2); end servicetonode(randoutnode(i),2)=b; servicetonode(randoutnode(i),3)=a; end end sumd=sum(servicetonode(

33、:,3); finalstnh=servicetonode; finalsd(h)=sumd;enda,v=min(finalsd)附錄5%a區(qū)交巡警服務(wù)平臺設(shè)置合理性指標(biāo)計算clear all;mln = xlsread(d:b,2,a2:b141); %matrix of line between nodesmn = xlsread(d:b,1,a2:e93); %matrix of nodes% create 0-1 matrix of line beteen nodesmln01=zeros(92,92);mlnsize=size(mln);mlnr=mlnsize(1); %row

34、number of mlnfor i=1:mlnr if mln(i,1)=92 & mln(i,2)30 path=path+(sp(fservices(j),i)-30)2; end end totalpath(i)=path; if work6.225 totalwork(i)=(work-6.225)2; largem=largem+1; endendsorttw,positionw=sort(totalwork)sorttp,positionp=sort(totalpath)附錄6%b區(qū)交巡警服務(wù)平臺設(shè)置合理性指標(biāo)計算clear all;load partb.mat % partb.

35、mat sp mn mln sn b131handlenode=minhn;len=length(mn);lensn=length(sn);avework=sum(mn(:,5)/lensn;totalwork=zeros(lensn,1);totalpath=zeros(lensn,1);largem=0;for i=1:lensn fservices=find(handlenode=sn(i); work=0; path=0; for j=1:length(fservices) work=work+mn(fservices(j),5); if sp(fservices(j),i)30 path=path+(sp(fservices(j),i)-30)2; end end totalpath(i)=path; if workavework totalwork(i)=(work-avework)2; largem=large